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models  were  numerically  integrated  through  time  by  a  linear  Euler  extrapolation 
technique.  A  varable  time  step  algorithm  was  included  that  maximized  time 
step  size  during  the  analysis  while  maintaining  good  accuracy.  This  program 
was  used  as  the  plane  stress  theoretical  model  for  the  HEN  procedure  to  analyze 
sustained  load  creep  crack  growth. 

A  method  for  getting  creep  crack  growth  behavior  solely  from  high  resolu¬ 
tion  displacement  measurements,  in  conjunction  with  a  cracked  specimen  model 
which  utilizes  realistic  constitutive  relationships,  has  been  developed. \JThe 
constitutive  law,  in  the  form  of  the  Bodner-Partom  material  model,  was 
cially  tailored  to  the  nickel-base  alloy  studied  which  displays  time  dependent 
nonlinear  inelastic  behavior  at  elevated  temperatures.  It  has  been  demonstrated 
that  the  technique  can  be  applied  where  crack  extension  is  very  small  and  co^uld 
not  otherwise  be  resolved  by  conventional  experimental  crack  measuring  tech-' 
niques.  This  method  provides  realistic  monotonically  increasing  crack  growth 
values.  The  predictions  agreed  to  within  10°>  of  post-test  measurements. 

Crack  growth  rate  and  crack  growth  criteria  were  studied.  Crack  tip 
strain  and  crack  opening  displacement  were  studied  in  the  HEN  results  for 
a  unique  parameter  controlling  crack  growth.  Because  the  parameters  were  not 
independent  of  time  due  to  apparent  environmental  degradation,  it  became 
necessary  to  establish  an  empirical  criterion  for  crack  growth  based  on  the 
best  fit  of  HEN  results.  A  damage  accumulation  criterion  based  on  creep 
rupture  formulations  was  also  developed  and  applied  with  promising  results. 

Several  crack  growth  rate  criteria  were  investigated,  one  of  which  is  the 
stress  intensity  factor.  The  K  criterion  matched  fairly  well  with  an  extrap¬ 
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SECTION  I 
INTRODUCTION 


1 .  BACKGROUND 

The  United  States  Air  Force  currently  places  strict  fracture 
mechanics  requirements  on  airframe  construction  and  maintenance  (Refer¬ 
ences  1,2).  These  requirements  involve  detection  of  flaws  by  periodic 
nondestructive  inspection  and  then  predicting  the  remaining  useful  life 
of  the  part  through  specified  fracture  mechanics  techniques.  Conse¬ 
quently,  if  an  airframe  part  is  examined  and  found  to  have  no  flaws 
that  can  grow  to  critical  size  prior  to  the  next  periodic  examination, 
it  can  be  returned  to  service. 

In  contrast  to  the  airframe,  low-cycl e- fatigue  limited  jet  engine 
parts  may  be  retired  from  service  when  no  flaws  have  yet  been  found  in 
them.  This  situation  occurs  because  the  retirement  of  engine  disks  is 
based  on  a  "crack  initiation"  criterion.  Under  this  criterion  all 
components  of  a  given  population  are  considered  to  have  failed  as  soon 
as  a  crack  of  some  finite  size  (e.g.,  .031  inches)  has  statistically 
formed  in  the  member  of  the  population  which  has  minimum  strength 
properties  (Reference  3).  No  attempt  is  made  to  utilize  the  additional 
life  associated  with  the  remaining  population  members  which  have 
statistically  higher  properties  and  are  expected  to  be  uncracked. 

From  a  safety  standpoint,  this  approach  has  been  generally  very 
successful.  But  for  real  materials  and  real  design  situations,  life¬ 
times  based  on  time  crack  initiation  of  the  minimum  member  tends  to  be 
extremely  conservative  for  a  component  population. 

It  has  been  estimated  that  replacement  costs  for  low-cycle- fatigue 
limited  jet  engine  disk  components  could  reach  the  $100,000,000  level 
by  the  1980  to  1985  time  period  (Reference  4).  A  significant  reduction 
of  this  cost  could  be  realized  if  a  procedure  was  developed  to  provide 
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accurate  assessment  of  useful  residual  life  in  retired  engine  disks. 

This  procedure  would  require  both  improved  inspection  and  fracture 
mechanics  life  prediction  techniques. 

The  present  research  is  aimed  at  developing  more  accurate  fracture 
mechanics  life  prediction  capabilities.  Various  controlling  aspects  of 
fatigue  crack  growth  are  strain,  stress,  stress  intensity,  temperature, 
load  application  frequency,  and  environment.  Speidel  (Reference  5) 
provides  an  excellent  discussion  of  the  relative  effects  of  each  of  the 
above  aspects  on  fatigue  crack  growth  rates  at  high  temperatures.  It 
was  shown  that  at  each  elevated  temperature  there  is  a  critical  frequency 
below  which  the  crack  growth  rate  is  creep  dependent  (i.e.  dependent  on 
exposure  time  to  load)  and  this  creep  dependency  increases  with  decreas¬ 
ing  frequency.  Also  the  effects  of  an  aggressive  environment  is  another 
time-dependent  phenomenon  that  results  in  a  frequency  dependence  of  the 
crack  growth  rate  which  is  very  similar  to  that  brought  about  by  creep. 

Various  engine  missions  may  include  long  dwell  times  at  high  stress 
levels.  Crack  growth  controlling  aspects  may  then  be  reduced  to  stress 
and  strain  levels,  stress  intensity,  temperature  and  environment  since 
load  cycling  is  not  occurring  during  these  dwell  periods. 

Due  to  the  high  stress  concentration  in  the  vicinity  of  the  crack 
tip,  (i.e.,  infinite  stress  concentration  using  linear  elastic  fracture 
mechanics)  and  the  high  temperature  environment  for  an  engine  disk,  the 
stress-strain  relations  for  the  material  are  nonlinear  and  time-dependent. 
The  high  stress  concentration  causes  the  material  to  yield  and  envelop 
the  crack  tip  with  a  plastic  zone.  Simultaneously,  the  high  temperature 
allows  the  material  under  load  to  flow  with  time,  the  phenomena  known 
as  creep,  (i.e.,  increase  strain  with  no  increase  in  stress).  Also  the 
environment,  both  temperature  and  atmosphere,  may  be  changing  the 
ductility  of  the  material  thereby  lowering  its  ability  to  strain  before 
fracture. 
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Creep  crack  growth  studies  of  metals  in  a  vacuum  indicate  that 
removal  of  air  from  the  crack  can  reduce  the  crack  growth  rates  by  two 
orders  of  magnitude  (Reference  5).  This  result  demonstrates  the 
importance  of  recognizing  environmental  effects  in  addition  to  what 
may  be  called  mechanical  effects  on  crack  growth.  Environmental  effects 
may  be  thought  to  dictate  the  critical  level  of  strain  at  fracture 
whereas  mechanical  effects  determine  the  rate  at  which  the  material 
moves  to  the  critical  strain  levels. 

The  present  work  mainly  deals  with  the  mechanical  aspects  of  crack 
growth  under  fixed  load  in  materials  that  deform  with  time  (creep  crack 
growth).  But  since  the  experimental  data  used  here,  came  from  elevated 
temperature  tests  in  laboratory  air,  some  environmental  effects  were 
also  considered. 


A  theoretical  model  for  creep  crack  growth  under  varying  and  fixed 
loads  must  be  able  to  account  for  the  changing  boundary  conditions 
associated  with  crack  growth.  In  the  event  of  total  unloading  and 
reloading  between  two  different  fixed  loads  the  possibility  of  crack 
closure  and  separation  needs  to  be  taken  into  account.  These  changing 
boundary  conditions  coupled  with  nonlinear  time  dependent  material 
behavior  are  well  suited  for  the  finite  element  method. 

2.  APPROACH 

The  modeling  effort  considered  here  involved  developing  a  two- 
dimensional  (plane  stress/plane  strain)  nonlinear,  time  dependent,  finite 
element  program  to  investigate  creep  crack  growth  under  constant  load. 

The  finite  element  analysis  incorporated  the  constant-strain  triangular 
element.  The  nonlinear  time-dependent  material  constitutive  model  took 
the  form  of  the  Bodner-Partom  viscoplastic  flow  law  (References  6, 7, 8, 9). 
This  flow  law  was  integrated  through  time  by  an  Euler  extrapolation 
scheme  (Reference  10)  and  incorporated  into  the  overall  finite  element 
program  by  means  of  the  residual  force  method  (Reference  11).  Material 
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constants  for  the  Bodner  material  model  were  determined  by  Stouffer 
(Reference  12)  to  best  match  the  behavior  of  IN-100  (Reference  13), 
a  current  jet  engine  turbine  disk  alloy. 

Time  step  sizes  of  the  Euler  scheme  were  maximized  subject  to 
specified  amounts  of  change  in  stress  and  strain  over  a  given  time  step. 
This  time  step  maximization  scheme  provided  the  ability  to  transition 
from  small  fraction  of  a  second  time  step  for  the  load  up  phase  to  large 
time  steps  of  the  order  of  minutes  or  even  hours  for  the  constant  load 
creep  phase.  This  variable  time  step  capability  was  a  necessity  to  make 
a  numerical  study  of  creep  crack  growth  computationally  feasible. 

Crack  growth  and  possible  crack  closure  during  unloading  was 
accounted  for  by  simple  modifications  to  the  structural  stiffness  matrix. 
These  simple  modifications  were  made  possible  by  choosing  an  iterative 
Gauss-Seidel  linear  equation  solver  (Reference  14)  which  requires  no 
explicit  factorization  of  the  stiffness  matrix.  Hence  between  time 
steps,  pertinent  terms  of  the  stiffness  matrix  could  be  easily  changed 
to  account  for  crack  growth  and  the  general  procedure  continued  without 
costly  matrix  factorization  time  required. 

The  finite  element  program  which  includes  the  capability  of  account¬ 
ing  for  material  creep  behavior  and  crack  growth  was  used  to  study  the 
creep  crack  growth  in  a  center  cracked  plate  test  specimen.  Several 
finite  element  models  were  incorporated  to  study  different  initial  crack 
lengths  in  the  plate  geometry.  These  models  were  subjected  to  various 
loads  that  were  chosen  to  coincide  with  a  parallel  experimental  program 
conducted  by  W.  Sharpe  (Reference  15). 

The  primary  objectives  in  the  present  research  were  to  determine 
the  actual  rate  of  creep  crack  growth  in  test  specimens  from  experimental 
displacement  and  compliance  measurements  and  to  determine  the  most 
reliable  criterion  for  predicting  creep  crack  growth  in  a  typical  jet 
engine  turbine  disk  alloy.  Specifically  the  interest  was  in  IN-100  at 
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1350°F.  These  objectives  were  to  be  accomplished  by  the  so-called  hybrid 
experimental-numerical  procedure  (Reference  16).  In  this  procedure, 
crack  growth  rates  would  be  estimated,  imposed  on  the  finite  element 
simulation,  and  then  the  finite  element  model  displacement  results 
would  be  compared  with  experimental  data  for  the  same  geometry  and 
loading  conditions.  A  similar  method  is  to  allow  the  crack  to  grow 
sufficiently  so  that  predicted  crack  opening  displacement  rates  from 
the  finite  element  model  match  experimental  data  for  the  same  geometry 
and  loading  conditions.  After  good  correlation  of  displacement.;  between 
model  and  experiment  was  achieved,  fracture  criteria  were  sought  out 
from  the  calculated  parameters  in  the  finite  element  model  such  as 
stress,  strain,  and  displacement.  Criteria  were  sought  which  could 
match  the  now  determined  experimental  crack  growth  rates,  displacement, 
or  displacement  rates  over  a  range  of  crack  length  and  load  levels. 

Once  a  reliable  creep  crack  growth  criterion  is  determined  and 
found  to  be  independent  of  specimen  geometry  it  can  then  be  applied  to 
an  actual  turbine  disk  specimen.  With  the  determination  of  flaw  sizes 
in  the  disk  through  nondestructive  examination  techniques,  the  crack 
growth  criterion  could  then  be  used  to  predict  the  remaining  time  for 
these  flaws  to  grow  to  critical  dimensions  under  constant  load 
applications. 
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SECTION  II 
LITERATURE  REVIEW 

Creep  crack  growth  in  general  may  be  thought  of  as  subcritical 
crack  growth  in  a  material  that  deforms  with  time  under  constant  external 
load.  This  time  dependent  deformation  or  creep  may  be  reversible 
(i.e.,  anelastic  creep)  or  it  may  be  permanent  (i.e.,  plastic  creep). 

At  elevated  temperature,  metals  generally  exhibit  nonlinear  time 
dependent  deformation.  Under  uniaxial  tensile  loading,  the  strain  in  a 
smooth  bar  increases  with  time  until  failure  ultimately  occurs.  Based 
on  similar  response  of  many  materials,  researchers  have  subdivided  the 
creep  curve  into  three  regions  as  shown  in  Figure  1.  After  the  initial 
instantaneous  strain  eQ,  materials  often  undergo  a  period  of  transient 
response  where  the  strain  rate,  e,  decreases  with  time  to  a  minimum 
steady-state  value  that  persists  for  a  substantial  portion  of  the 
materials  life.  These  two  regions  are  referred  to  as  transient  or 
primary  stage  and  steady-state  or  secondary  stage  respectively.  Final 
failure  with  rupture  life  tR  occurs  soon  after  the  creep  rate  begins 
to  increase  during  the  third  or  tertiary  stage  of  creep.  A  common 
empirical  relationship  between  creep  strain  rate  and  stress  in  the 
secondary  stage  of  creep  is  given  as: 

e  =  yc  (a)S  (1) 

where  a  is  the  uniaxial  stress,  e  is  the  creep  strain  rate,  and  y 

c 

and  6  are  empirical  constants  chosen  to  match  creep  test  results. 

Creep  crack  growth  has  been  studied  using  viscoelastic  (References 
17-20),  viscoelastic-plastic  (Reference  21),  and  plastic  creep  material 
models  for  metals  as  indicated  in  a  literature  review  by  Fu  (Reference 
22).  The  viscoelastic  modeling,  a  form  of  anelastic  behavior,  pertains 
mainly  to  nonmetals  such  as  elastomers,  polymers,  and  solid  rocket 
propellants.  The  material  of  interest  in  this  investigation  is  a 
current  jet  engine  turbine  disk  alloy  known  as  IN-100  (Reference  13). 
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Figure  1.  Three  Stages  of  Creep  Behavior 


Figure  2.  Line  Integra]  Paths  for  J  and  C* 
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In  the  present  study  all  time  dependent  deformation  of  IN-100  is  treated 
as  permanent  and  consequently  any  anelastic  behavior  is  considered  negli¬ 
gible.  However,  a  brief  review  of  some  anelastic  (i.e.  viscoelastic) 
creep  crack  growth  studies  will  be  addressed  later. 

The  next  two  review  sections  will  consider  recent  creep  crack  growth 
studies  of  an  analytical  or  closed  form  nature  as  well  as  applications 
of  the  finite  element  method  to  creep  crack  growth. 

1.  ANALYTICAL  CREEP  CRACK  GROWTH  STUDIES 

Knauss  (Reference  17)  analytically  modeled  steady  crack  growth  in  a 
viscoelastic  sheet.  In  his  study  the  plastic  zone,  which  was  assumed 
small,  was  accounted  for  by  prescribed  fixed  and  finite  stress  distri¬ 
butions  in  the  crack  tip  region.  No  interaction  between  crack  tip  and 
the  far  field  stresses  were  allowed.  This  means  the  stress  profile 
ahead  of  the  moving  crack  tip  remained  constant  and  independent  of  the 
far  field  stresses.  Magnitudes  and  gradients  of  the  stress  in  this 
crack  tip  region  were  studied  along  with  two  crack  growth  criteria. 

The  two  criteria  were  maximum  strain  and  a  maximum  strain  energy 
cri terion. 

Schapery  (References  18,19,20)  performed  a  viscoelastic  crack  growth 
analysis  similar  to  Knauss  but  placed  no  significant  restrictions  on 
the  nature  of  the  failing  material  at  the  crack  tip.  It  could  be  highly 
nonlinear  and  rate  sensitive.  An  energy  criterion  for  failure  was  also 
used  here. 

Wnuk  (Reference  21)  included  plasticity  with  viscoelasticity  for 
his  quasistatic  extension  of  a  tensile  crack  analysis.  A  "final  stretch" 
crack  growth  criterion  was  proposed.  This  criterion  postulates  that  the 
amount  of  deformation  which  occurs  within  the  crack  tip  region  or  pro¬ 
cess  zone  during  the  time  interval  just  prior  to  decohesion  of  this  zone 
is  a  material  constant.  In  contrast  to  the  maximum  strain  criterion, 
the  "final  stretch"  criterion  is  path-dependent. 
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Fu's  review  (Reference  22)  of  quasi-static  crack  growth  in  metals 
at  elevated  temperature  presents  four  different  creep  crack  growth  rate 
equations: 


a 


a  =  a  (k) 

(2) 

a  =  B  (cn)0 

(3) 

a  =  C  (y)n 

(4) 

a  =  D  (C*)* 

(5) 

where  a  is  the  crack  growth  rate,  K  is  the  linear  elastic  stress  inten¬ 
sity  factor,  on  is  the  net  section  stress,  (e.g.,  load  divided  by  remain¬ 
ing  uncracked  ligament  in  a  center  cracked  plate  geometry),  y  in  general 
is  the  load-point  displacement  rate,  and  C*  is  a  line  integral  related 
to  the  rate  of  change  of  potential  energy  release  per  unit  of  crack 
growth  (References  22-25).  The  C*  integral  also  discussed  in  Appendix 
C,  is  obtained  directly  from  Rice's  J  integral  by  introducing  strain 
rate  and  displacement  rate  instead  of  strain  and  displacement  such  that: 


y.  9U  . 

0  =  f  Whjto  -  T,  ds] 


becomes 


j.  3U . 

C*  =  f  [W*(ei;j)dy  -  Ti  ^  ds] 


where 


/mn 

a .  .  de  .  . 
1J  1J 


“*('u)  f  m  ”ij  dEij 


7..  =  traction  vector 


u..,  u^  =  displacement  and  displacement  rate  respectively 
e,.  =  strain  and  strain  rate  respectively 

1  J  '  J 


(6) 


(7) 


(8) 

(9) 


r,  x,  ds  -  see  Figure  2 
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There  are  available  experimental  data  supporting  any  of  the  rate,  see 
Equations  2  through  5  as  listed  by  Fu. 

The  following  conclusions  have  been  extracted  from  the  literature: 

1.  The  stress  exponent  3  in  Equation  1  plays  an  important  role  in 
determining  the  characterization  of  the  crack  tip  behavior.  For 

3  <_  5  the  stress  intensity  factor  approach,  Equation  2  may  be  used, 
and  for  6  >  7  the  net  section  stress  approach  may  be  used  (Refer¬ 
ences  26,27). 

2.  Critical  test  conditions  for  evaluating  a  creep  crack  growth 
criterion  should  consist  of  at  least  two  geometries  which  have 
different  stress-intensity-factor  divided-by-stress  ratios  (Refer¬ 
ence  23).  Some  crack  growth  criteria  have  been  found  to  be  depen¬ 
dent  on  geometry  and  therefore  have  no  general  application.  Use 
of  two  or  more  test  geometries  helps  determine  how  dependent  a 
crack  growth  criterion  is  on  geometry. 

3.  Creep  crack  growth  results  from  two  competing  processes.  These 
processes  are:  (1)  growth  and  coalescence  of  defects  which  contri¬ 
bute  to  crack  advancement  and  (2)  the  creep  deformation  process 
that  causes  retardation  and  even  arrest  of  crack  growth  (Reference 
25). 

4.  Creep  crack  growth  rates  are  very  sensitive  to  environmental 
effects.  Removal  of  air  from  the  crack  can  reduce  the  crack  growth 
rates  by  two  orders  of  magnitude  (Reference  5). 

5.  Crack  opening  displacement  crack  growth  theories  indicate 
that  failure  times  due  to  creep  crack  growth  are  controlled  by 
the  stress  intensity  factor  at  large  stresses  and  by  net  section 
stresses  at  very  low  stresses  (Reference  28).  However,  a  counter 
viewpoint  is  stated  in  Reference  29  where  it  is  concluded  that 
creep  crack  growth  does  not  correlate  well  with  the  stress  intensity 
factor  at  relatively  high  stress  levels. 
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6.  Creep  crack  growth  rates  correlate  with  the  energy  rate  integral 
C*.  This  method  holds  great  promise  for  design  calculations  because 
C*  can  be  calculated  using  finite  element  analysis  as  well  as 
measured  empirically  in  constant  displacement  rate  tests  (Reference 
24). 

7.  An  approximation  to  C*  in  the  form  of  the  product  of  net 
section  stress  and  load  line  displacement  rate,  referred  to  as  3, 
gives  good  correlation  of  creep  crack  growth  rates  in  specimens 
of  different  geometries  (Reference  29). 

8.  Crack  growth  theories  generally  fall  into  one  of  two  categories. 
Either  they  are  of  an  energy  nature  (e.g.,  J  or  C*  integrals),  or 
they  deal  with  some  localized  crack  tip  parameter  such  as  strain 

or  crack  opening  displacement. 

2.  FINITE  ELEMENT  ANALYSIS  RELATED  TO  QUASI-STATIC  TIME  DEPENDENT 
CRACK  GROWTH 

The  general  technique  of  approximating  a  continuum  with  simple 
discrete  elements  such  as  uniaxial  bar  elements  dates  back  to  the  1940's 
(Reference  30).  The  history  of  the  finite  element  method  in  structural 
analysis,  as  it  is  known  today,  is  well  described  by  Zienkiewicz  (Refer¬ 
ence  31).  The  application  of  the  method  to  the  problem  of  nonlinear 
material  behavior  has  also  been  developed  (Reference  31).  The  purpose 
of  this  section  is  to  briefly  review  the  use  of  the  finite  element 
method  for  the  stress  analysis  of  cracked  plates  where  nonlinear  time 
independent  and  time  dependent  materials  models  were  employed. 

a.  Elastic-Plastic  (Time  Independent)  Analysis 

The  finite  element  method  has  been  widely  used  to  determine  the 
stress  and  strain  fields  around  cracks  in  nonlinear  materials  where 
time  independent  elastic-plastic  materials  models  were  incorporated 
(References  32-45).  Several  of  these  use  the  "initial  stress"  or 
"initial  strain"  approach  t.o  elastic-plastic  modeling  as  described  in 
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Appendix  A.  During  the  last  few  years  many  reports  have  been  published 
on  this  subject.  This  review  will  not  cover  all  these  publications, 
but,  it  will  concentrate  on  two  areas:  (1)  crack  tip  element  selection 
and  (2)  crack  growth  modeling  by  the  finite  element  method. 

There  are  several  choices  to  consider  when  selecting  finite  elements 
to  model  the  vicinity  of  a  crack  tip.  These  choices  might  best  be 
classified  into  two  categories.  The  first  category  is  to  use  the  same 
element  type  incorporated  in  the  model  remote  from  the  crack  tip  (e.g., 
constant  strain  triangles).  The  second  category  is  to  use  a  special 
crack  tip  element  that  has  functions  to  include  the  specific  stress  or 
strain  singularity  desired  at  the  crack  tip. 

There  are  of  course  benefits  and  disadvantages  to  each  choice. 

Using  the  same  element  such  as  a  constant  strain  triangle  both  at  the 
crack  tip  and  remote  definitely  has  the  benefit  of  simplicity.  However, 
the  cost  of  this  simplicity  is  the  requirement  to  use  large  numbers  of 
such  elements  around  the  crack  tip  to  obtain  acceptable  results.  Also 
in  a  fully  plastic  material,  some  elements  do  not  accurately  model 
incompressible  strain  behavior  (Reference  46).  The  bilinear  displace¬ 
ment  quadralateral  is  most  susceptible  to  this  inaccuracy,  furthermore 
the  constant  strain  triangle  is  one  of  the  least  susceptible. 

The  selection  of  a  special  crack  tip  element  requires  knowledge 
of  the  stress  or  strain  singularity  around  the  crack  tip.  Several 
authors  have  developed  the  inverse  square  root  "r"  singularity  for 
linear  elastic  crack  tip  behavior.  Accordingly  many  special  elements 
have  been  addressed  to  this  singularity.  In  a  paper  entitled  "Crack 
Tip  Finite  Elements  Are  Unnecessary"  (Reference  47),  the  authors 
describe  the  modification  of  the  eight  noded  isoparametric  element  such 
that  it  incorporated  the  1//T  stress  singularity.  Therefore,  any 
existing  finite  element  program  tnat  had  the  eight  noded  isoparametric 
element  in  it  also  effectively  has  one  form  of  a  special  crack  tip 
element  capability  for  elastic  analysis. 
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where  K  ,  K  are  plastic-intensity  factors  and  n  is  a  strain  hardening 
exponent.  This  type  of  singularity  was  embedded  in  a  crack  tip  element 
and  used  to  evaluate  the  plastic  intensity  factors  as  a  function  of 
applied  loading  (Reference  49).  The  advantage  of  this  method  was  that 
the  large  strain  gradients  in  the  crack  tip  region  are  accounted  for  by 
Equation  10.  For  elastic-plastic  and  creep  type  materials  the  singu¬ 
larity  changes  between  initiation  and  growth  of  the  crack  and  in  general 
is  not  known.  The  development  of  a  special  crack  tip  element  in  the 
case  where  the  singularity  is  not  known  to  begin  with  would  be  an 
extensive  undertaking  by  itself  (Reference  37)  without  bringing  in  the 
additional  complexity  of  a  moving  crack  tip. 


Finite  element  researchers  have  considered  the  crack  growth  problem. 
Kobayashi ,  Chiu,  and  Beeukes  (Reference  43)  analyzed  an  extending  crack 
under  monotonically  increasing  load.  Crack  extension  was  achieved  by 
applying  a  relief  force  equal  in  magnitude  but  opposite  in  direction 
to  the  restraining  force  at  the  crack  tip  node.  This  relief  force  was 
applied  in  100  equal  increments  or  in  one  single  increment.  The  crack 
opening  displacements  at  the  node  adjacent  to  the  crack  tip  computed 
by  the  single  increment  method  were  less  than  5  percent  smaller  than 
the  corresponding  displacements  computed  with  the  100  increment  method. 
Lee  and  Liebowitz  (Reference  35)  using  a  similar  technique  demonstrated 
that  plastic  strain  energy  increased  linearly  with  crack  length  as  the 
crack  grew  in  their  analysis.  Anderson  (Reference  36)  also  made  use  of 
relaxing  the  crack  tip  node  force  incrementally  to  simulate  crack  growth 
in  the  finite  element  model. 
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Shi h ,  DeLorenzi,  and  Andrews  (Reference  42)  analyzed  crack  initia¬ 
tion  and  stable  crack  growth  in  elastic-plastic  material  by  special  use 
of  eight-node  i soparametric  elements.  Initial  crack  tip  blunting  was 
modeled  followed  by  crack  extension.  Crack  extension  was  modeled  by  a 
combination  of  crack  tip  node  shifting  and  then  releasing  it  to  move 
on  to  the  next  element's  nodes.  Moving  these  nodes  within  elements  at 
the  crack  tip  required  approximations  to  be  made  about  Gauss  point 
location  and  relocation.  These  approximations  are  not  necessary  when 
only  a  node  release  method  and  elements  such  as  constant  strain  triangles 
are  employed  to  model  crack  growth.  The  1/r  strain  singularity  provided 
by  the  special  use  of  these  eight  noded  isoparametric  elements  may  better 
model  the  crack  tip  singularity  for  a  fixed  crack  length  and  a  theoreti¬ 
cal  continuum.  However,  when  the  crack  initiates  and  grows  in  a  creep 
type  material  the  strain  singularity  is  unknown,  especially  when  consider¬ 
ing  a  grain  structure  around  the  crack  tip  rather  than  a  continuum  and 
the  fact  that  a  creep  crack  follows  an  intergranular  path. 

b.  Visco-Plastic  and  Creep  (Time  Dependent)  Analysis 

Several  references  have  been  found  on  viscous  or  time  dependent 
material  models  being  incorporated  into  finite  element  programs  (Refer¬ 
ences  50-60).  In  general  these  material  models  may  be  similar  to 
Equation  1  for  pure  creep  with  the  addition  of  time  independent  elastic- 
plastic  relationships,  or  they  may  have  sho^t-term  response  viscoplastic 
relationships  that  only  model  the  load  up  phase.  Zienkiewicz  (Reference 
10)  proposed  placing  a  short-term  response  viscoplastic  equation  in 
series  with  a  long-term  response  creep  law,  similar  to  Equation  1. 

This  would  be  a  unified  time  dependent  material  model  where  no  direct 
coupling  is  assumed  between  short-term  plastic  strains  and  long-term 
creep  strains. 

Only  a  few  papers  have  been  found  to  date  that  use  a  time  dependent 
material  model  and  the  finite  element  technique  to  model  crack  growth. 
Ohtani  ano  NaKamara  (Reference  61)  analyzed  creep- crack,  propagation  with 
an  elastic-secondary  creep  material  model.  The  secondary  creep  law  was 


AFWAL-TR-80-41 40 


identical  to  Equation  1  in  its  uniaxial  form.  A  critical  strain  criterion 
(i.e.  average  effective  plastic  strain  at  the  crack  tip)  was  used  for 
determining  the  time  to  grow  the  crack  in  the  finite  element  model  of  a 
center  cracked  plate.  Crack  tip  nodes  were  released  to  simulate  crack 
growth.  Goodall  and  Chubb  did  a  similar  analysis  on  the  compact  tension 
specimen  (Reference  62).  A  critical  strain  crack  growth  criterion  was 
again  used.  Finally  Zaphir  and  Bodner  (Reference  63)  incorporated  the 
time  dependent  Bodner-Partom  viscoplastic  flow  law  into  the  NONSAP  finite 
element  code  to  analyze  a  double-centi lever-beam  cracked  geometry.  In 
this  case,  high  loading  rates  were  studied  over  short  time  periods  rela¬ 
tive  to  creep  analyses.  Consequently  the  "recovery  term",  as  described 
in  a  later  section  and  used  to  best  model  creep,  was  not  included.  No 
crack  growth  was  allowed  here. 

In  each  of  the  finite  element  solutions  referred  to  previously, 
a  unified  time-dependent  material  model  that  not  only  accurately  models 
the  short-tern:  load  up  stage  of  material  response,  but  also  transition 
smoothly  into  the  pure  creep  stage  was  never  considered.  In  addition, 
the  hybrid  experimental-numerical  technique  was  not  used  with  high 
resolution  experimental  crack  opening  displacement  data.  It  was 
anticipated,  in  the  present  research,  that  the  combination  of  a  more 
realistic  time  dependent  material  model  and  high  resolution  test  data 
would  result  in  a  much  better  understanding  of  what  controls  creep 
crack  growth  than  provided  by  these  prior  analyses. 
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SECTION  III 

FINITE  ELEMENT  COMPUTER  PROGRAM  DEVELOPMENT 

This  section  describes  the  development  and  validation  of  a  finite 
element  computer  program  called  "VISCO".  VISCO  is  a  two-dimensional 
plane  stress/plane  strain  program  that  incorporates  three  different 
nonlinear  time  dependent  viscoplastic  material  models.  It  uses  constant 
strain  triangular  elements  with  the  option  to  release  fixed  nodes  and 
thus  has  the  capability  to  simulate  crack  growth.  Results  from  VISCO 
are  compared  with  other  published  solutions  to  check  its  validity. 

The  approach  selected  for  elastic-viscoplastic  analysis  with  the 
finite  element  technique  employs  the  "residual  force  method"  (Reference 
11)  (see  Appendix  A  for  a  complete  discussion  on  this  particular 
approach).  In  the  residual  force  method  the  elastic  stiffness  is  used 
during  the  entire  analysis  and  any  nonlinear  elastic-viscoplastic  deforma¬ 
tion  that  occurs  must  be  accounted  for  by  developing  so-called  plastic¬ 
load  vectors  that  add  to  the  force  side  of  the  governing  equilibrium 
equation.  In  general,  the  matrix  equation  which  governs  the  response 
of  a  discretized  structure  can  be  written  as 

[K]  {U}1  =  {P}1  +  {Q}1'"1  (11  ) 

where  [K]  is  the  elastic  stiffness  matrix  {U}1  is  the  generalized  nodal 
displacement  vector  for  the  ith  time  step,  {P}1  is  the  load  vector  after 
the  ith  time  step  due  to  external  forces,  and  {Q}1"1  is  the  plastic-load 
vector  computed  from  plastic  strains  accumulated  prior  to  the  ith  time 
step.  For  each  element,  these  plastic-load  vectors  are  self-equilibrat¬ 
ing.  The  viscoplastic  strain  rate  expressions  which  develop  plastic 
strains  with  time  under  load  are  described. 

1.  VISCOPLASTIC  MATERIALS  MODELS 

In  solid  mechanics  it  is  customary  to  separate  the  two  important 
groups  of  phenomena  described  respectively  by  "creep"  and  "plasticity". 


is 
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The  first  includes  all  time  dependent  effects  and  results  in  creep 
strains  accumulating  at  a  finite  rate.  The  second  group  develops 
permanent  (plastic)  strains  Instantaneously  since  time  does  not  enter 
directly  into  consideration  as  in  the  elastic-plastic  approaches  given 
in  Appendix  A.  From  a  physical  point  of  view,  creep  and  plasticity 
cannot  be  treated  separately  as  only  the  combined  effect  is  measurable. 
Also  the  concept  of  time  independent  or  instantaneous  plasticity  is 
only  a  convenient  mathematical  approximation  and  not  experimentally 
based. 


Viscoplasticity,  in  a  complete  sense,  is  the  combination  of  these 
two  strain  groups  into  a  unified  plastic  strain  rate  model.  A  model 
with  this  capability  is  the  Bodner-Partom  viscoplastic  flow  law  (Refer¬ 
ences  6-9)  from  hereon  referred  to  as  the  Bodner  model.  The  superposi¬ 
tion  of  Malvern's  overstress  law  (References  66,67)  with  Norton's  law 
for  secondary  creep  has  also  been  proposed  as  a  unified  viscoplastic 
flow  law  (Reference  10). 


Each  of  these  flow  laws  has  been  incorporated  into  VISCO  by  assum¬ 
ing  small  strains  and  decomposing  the  total  strain  rate  into  elastic 
(reversible)  and  plastic  (nonreversibl e)  components. 


•  _  -E  ,  -P 

e  .  ■  -  e .  .  +  e  .  . 
1 J  1J  1J 


(12) 


which  in  general  are  both  nonzero  for  all  loading/unloading  conditions. 
Anelastic  stresses  and  strains  corresponding  to  time  dependent  reversible 
deformations  with  energy  losses  are  not  considered  in  this  formulation 
and  are  assumed  to  be  relatively  unimportant.  This  assumption  seems 
quite  justified  based  on  the  good  correlation  between  predicted  and  test 
results  to  be  shown  later. 
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The  elastic  strain  rate  e..E  is  related  to  the  stress  rate  by  the 

*  J 

p 

time  derivative  of  Hooke's  law.  The  plastic  strain  ,  assuming 

Incompressibility  and  isotropy,  is  taken  to  follow  the  Prandtl-Reuss 
flow  law  of  classical  plasticity 


(13) 


where  are  the  components  of  the  deviatoric  stress  tensor  and  A  is  a 
scalar  that  reflects  the  viscosity  of  the  material.  The  specific  form 
of  A  is  presented  below  for  each  of  the  viscoplastic  flow  laws. 


a.  Malvern  (Overstress)  Flow  Law 


A  portion  of  a  total  viscoplastic  model  that  accounts  for  the  so- 
called  instantaneous  plasticity  during  loading  might  take  the  form  as 
given  by  Malvern  (References  66,67),  otherwise  known  as  the  "overstress" 
model : 
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(14) 


where  Yp  is  a  fluidity  constant  whose  magnitude  will  determine  the  strain 

rate  sensitivity  of  the  model,  see  Figure  3,  og  is  the  effective  stress 

defined  as  3  J2  where  J2  is  the  second  invariant  of  the  deviatoric  stress 

defined  as  J2  =  l/2S^jS.j,  and  o  (e£)  is  the  strain  hardening  yield  stress 

shown  to  be  a  function  of  the  effective  plastic  strain,  e£,  defined 

incrementally  as  de*3  = de.?  de?..  The  strain  hardening  stress  function, 

e  1  3  1J  1J 

within  VISCO,  takes  one  of  two  forms,  either  a  linear  relationship  such  as 

oUj)  -  a0  +  H'  *  EP  (15) 

or  a  Ramberg-Osgood  type 
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STRAIN 


Figure  3.  Variation  of  Malvern  Model  Behavior  with  Yp  Under  Constant 
Strain  Rate  Loading  Conditions 
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R(eP)  m  if  R(ep)m  >  o 

6  6  i  "  o  (16) 

o  if  R(eP)m  <  o 

o  e  o 

where  aQ  1S  the  initial  yield  stress,  the  constant  H‘  represents  the 
slope  of  the  stress  versus  plastic  strain  curve,  and  R  and  m  are  con¬ 
stants  of  a  Ramberg-Osgood  type  stress-strain  curve.  In  Equation  16, 

note  that  if  the  effective  plastic  strain  ep  is  at  or  near  zero  such 

e 

p  -m  _ 

that  the  function  R(e£)  is  less  than  the  initial  yield  stress,  oq, 

then  o(ep)  is  set  equal  to  o  . 

Implementation  of  the  Malvern  model  then  requires  selection  of 
Equation  15  or  Equation  16.  This  selection  would  depend  on  the  best  fit 
of  the  material's  uniaxial  stress-strain  Curve  developed  under  strain 
rates  at  or  near  the  lowest  strain  rate  expected  to  be  modeled  with  the 
Malvern  Law.  If  Equation  15  were  chosen,  the  initial  yield  stress,  aQ, 
and  slope,  H1  would  be  determined  from  this  experimental  curve.  A 
similar  determination  would  be  done  if  Equation  16  were  selected.  The 
fluidity  constant,  y  ,  would  be  chosen  to  best  reflect  the  strain  rate 
sensitivity  of  the  material  (see  Figure  3)  displayed  by  experimental 
stress-strain  data  developed  at  high  strain  rates. 

The  Malvern  model  may  also  be  used  to  perform  time  independent 
elastic-plastic  solutions.  In  this  case  time  becomes  a  fictitious 
parameter  and  thus  the  fluidity  constant,  y  ,  may  take  on  any  nonzero 
positive  value.  The  elastic-plastic  solution  is  che  steady  state  value 
of  the  stresses,  strains,  and  displacements  after  the  load  has  been 
applied.  This  has  been  found  to  occur  in  approximately  30  time-steps 
after  maximum  load  application  unless  total  section  yielding  develops. 
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An  Euler  linear  extrapolation  scheme  is  employed  within  VISCO  for 
the  time  integration  of  the  viscoplastic  strain  rate  expressions.  The 
Malvern  model  is  integrated  in  VISCO  as  follows 


•  P 

{eij}  = 


->] 


3  {!ii> 


i-1 


if  oe1_1  <  1 


aei_1  if  aei_1  >  a^ee^ 


i-1  (17) 


{de?.}1  =  {e?.}1  dt1 

■  J  *  J 


(18) 
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or  •  1 
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where  the  superscript  i  refers  to  the  time-step  and  a  subscript  i  refers 
to  specific  components  of  stress  or  strain. 


b.  Norton's  Law  for  Secondary  Creep 


Another  portion  of  a  total  viscoplastic  flow  law  that  complements 
the  Malvern  model  and  accounts  for  long-term  creep  is  given  by  Norton's 
creep  law  (Reference  68)  written  in  multiaxial  form  as 


•  P  _  ,  3  $.. 

eij  -  Yc  °e  2 

e 


(21) 


where  yc  and  8  are  constants  determined  from  uniaxial  creep  test  results. 

Creep  test  data  at  two  different  stress  levels  are  required.  A  straight 
line  is  fitted  to  each  test's  secondary  stage  of  creep  strain  plotted 
versus  time  (Figure  1).  The  slope  of  this  line  or  strain  rate  and  stress 
level  from  each  test  is  substituted  into  Equation  21  which  provides  two 
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equations  to  solve  for  the  two  unknown  constants.  After  taking  the 
natural  logarithm  of  Equation  21,  the  simultaneous  solution  of  the  two 
equations  is  straightforward. 


The  Euler  extrapolation  scheme  for  integrating  Equation  21  in  VISCO 
is  similar  to  the  one  given  for  the  Malvern  model  and  will  not  be 
repeated  here. 

It  appears  that  the  combination  of  the  Malvern  and  Norton  laws 
would  provide  a  complete  viscoplastic  flow  law.  But,  metals  at  elevated 
temperature  have  been  observed  to  display  a  phenomenon  called  "recovery" 
(Reference  69).  Recovery  is  the  softening  of  cold-worked  metal  or  it 
may  be  characterized  as  a  fading  memory  of  prior  strain  hardening. 


In  creep  crack  growth  at  elevated  temperature  and  under  constant 
load,  consider  material  well  ahead  of,  but  in  the  path  of  the  crack. 
During  load  up  this  material  will  plastically  deform  and  strain  harden 
depending  on  its  proximity  to  the  initial  crack  tip.  Later,  during 
the  sustained  load  phase,  the  phenomena  of  recovery  will  allow  this 
material  to  soften  prior  to  the  arrival  of  the  crack  tip  at  which  time 
additional  loading  will  occur.  The  amount  of  recovery  prior  to  the 
arrival  of  the  crack  tip  should  then  have  some  effect  on  the  values 
of  stress  and  strain  developed  around  the  crack  tip  when  it  reaches 
the  material  being  considered.  Therefore,  this  investigation  will 
focus  on  the  following  viscoplastic  flow  law  which  does  include  the 
phenomenon  of  recovery. 


c.  Bodner  Viscoplastic  Flow  Law 

In  this  formulation  by  Bodner  and  Partom  (References  6-9)  the  A 
parameter  from  Equation  13  is  expressed  in  terms  of  second  invariants 
by  making  use  of  the  square  of  Equation  3 


1  •  P 

2  e1j 


—  A^S  S 

2  A  '  J  i  j 


(22) 


p 

where  Dp  is  the  second  invariant  of  the  plastic  strain  rate  and  Jp  is 
the  second  invariant  of  the  deviatoric  stress. 
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Rather  than  specifying  a  specific  yield  criterion  as  in  classical 
plasticity,  this  formulation  by  Bodner  and  Partom  is  based  on  the  assump¬ 
tion  that  a  continuous  functional  relationship  exists  between  the  plastic 
deformation  rate  and  the  stress  invariants,  i.e. 

°2  =  D2(J2»Zk’T)  (23) 

where  Z^  are  one  or  more  internal  (viscoplastic)  state  variables  and  T 
is  the  temperature.  Introducing  Equation  23  into  Equation  22  and  solving 
for  x  gives 

P  1 

X  =  [D^J^Z^T)/  J2]  2  (24) 

The  general  form  for  the  evolution  equation,  i.e.  history  depen¬ 
dence,  of  the  viscoplastic  state  variables  is 

Zk  -  Fk(J2.Zk,T)  (25) 

For  conditions  cf  uniaxial  stress  of  constant  sign,  the  hardened 
state  with  respect  to  plastic  flow  is  assumed  to  be  represented  by  a 
single  state  variable  Z  which  depends  on  plastic  work.  This  single  state 
variable  Z  also  corresponds  to  isotropic  hardening.  Additional  state 
variables  are  necessary  for  such  characteristics  as  kinematic  hardening 
which  will  not  be  employed  here. 

D 

The  particular  form  cnosen  for  D-,  (^.Z.T)  was  motivated  by  the 
equations  of  dislocation  dynamics  and  given  by  Bodner  and  Partom  as 

3?  =  °o  'r-l  <26> 

The  factor  (n  +  1  )/n  was  introduced  at  an  early  stage  in  the  development 
of  the  equations  for  numerical  purposes  and  only  affects  the  numerical 
values  of  some  of  the  material  constants.  The  constant  is  described 
as  the  limiting  value  of  the  plastic  strain  rate  in  shear.  Its  value  can 
be  arbitrarily  chosen  and  is  usually  taker,  to  be  the  same  large  number 


23 


S 


AFWAL-TR-80-41 40 


4  -1 

for  all  materials.  A  value  DQ  =  10  sec  is  generally  adequate  except 
for  conditions  of  very  high  rates  of  straining. 

The  parameter  n  controls  strain  rate  sensitivity  and  also  influences 
the  overall  inelastic  level  of  the  stress-strain  curves.  A  decrease 
of  n  leads  to  increasing  strain  rate  sensitivity  and  lowering  of  the 
level  of  the  stress-strain  curves. 


As  stated  previously,  Z  is  assumed  to  be  a  function  of  plastic  work, 
Wp,  and  the  following  relationship  is  introduced 

Z  =  Z  (Wp)  =  Z1  -  (Z1-ZQ)  exp[-m  Wp]  (27) 

The  quantity  Z-|  is  the  maximum  value  of  Z  which  is  necessary  if  the 

deformations  do  not  revert  to  elastic  behavior  at  large  values  of  Wp. 

Z  is  the  value  of  Z  for  which  W„  =  o  and  can  therefore  be  the  initial 
o  p 

state  point  from  which  plastic  work  is  measured.  It  is  noted  that  the 
general  function  (Equation  27)  would  be  a  basic  material  property  and 
that  Wp  is  the  relative  amount  of  plastic  work  done  from  some  initial 

state,  (i.e.,  Wp  is  not  an  absolute  parameter). 

“p  a/sij  W dt  (28> 

The  quantity  m  in  Equation  27  is  a  material  constant  that  relates  to 
the  rate  of  work  hardening. 


At  high  temperatures,  it  is  generally  necessary  to  consider  the 
thermal  recovery  of  hardening  generated  by  plastic  deformation.  In 
this  case  the  plastic  work,  Wp,  is  redefined  as  follows: 


W 

P 


/Wo 


rec 

radT"-  Z) 


(29) 
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where 


(30) 


where  Z^  is  the  stable,  (i.e.,  non-work  hardened)  value  of  Z  at  a  given 

temperature,  and  A  and  r  are  additional  material  constants  chosen  to 
match  the  models  behavior  to  creep  test  data  as  was  done  for  the  Norton 
model.  Note  that  the  second  term  on  the  right  of  Equation  29  (i.e.,  the 
recovery  term)  makes  a  negative  contribution  to  Wp  due  to  the  negative 

sign  on  A,  since  Z  is  always  greater  than  or  equal  to  Z^. 


The  recovery  term  in  Equation  29  is  essential  if  the  material  model 
is  to  be  able  to  represent  secondary  creep.  Secondary  creep  is  the 
balanced  condition  when  the  rate  of  work  hardening  equals  the  rate  of 
thermal  recovery  or,  setting  the  time  derivative  of  Equation  29  to  zero. 


S.  .e.  .  + 

ij  iJ 


rec 

m(ZrZ) 


=  0 


(31) 


At  relatively  high  strain  rates,  the  thermal  recovery  in  Equation 
29  is  relatively  unimportant  and  the  steady  state  condition  is  realized 
when  Z  reaches  its  saturation  value  ly 


Again,  VISCO  employs  the  Euler  extrapolation  scheme  for  the  numerical 
time  integration  of  the  Bodner  equations.  During  each  time  step  Equation 
13  and  Equations  24  through  30  are  integrated  as  follows  for  each  element 

Z1  =  Z]  -  (ZrZ0)  exp[-m  Wp1_1]  (32) 

(02P)j  =  D02  exp[-  (  3^)"  (33) 

=  [(D^)1/^1"1]^  {S.j}1'1  (34) 

Ide.J}1'  =  (e.j  j)1  dt1  (35) 
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•  i 
Z 


rec 


W  =  W 
P  P 


i-1 


* 1 V 


i-1 


{de;.}  +  z  a 

i  j  rec 


dt1/[m(Z1-Z1)] 


(36) 

(37) 


where  the  superscript  i  refers  to  the  time-step  and  a  subscript  i  refers 
to  a  specific  material  constant,  Z.. ,  or  to  specific  components  of  stress 
and  strain. 


The  specific  material  constants  required  for  the  Bodner  model  in  this 
investigation  were  determined  to  best  fit  the  behavior  of  IN-100  at  the 
temperature  of  1350°F.  This  material  was  characterized  by  performing 
uniaxial  tensile  stress-strain  tests  and  creep  tests  at  different  stress 
levels  at  1350°F.  The  following  constants  for  the  Bodner  model  were 
developed  by  Stouffer  (Reference  12)  and  the  details  of  this  procedure 
are  summarized  in  Appendix  D.  The  A  and  r  constants  defined  by  Stouffer 
are  different  than  the  constants  used  in  the  present  formulation: 

D„  =  104  sec"1 
o 

n  =  3.50 
ZQ  =  224.4  ksi 
Z]  =  251.5  ksi 
m  =  3.750  ksi'1 

This  first  group  of  constants  is  based  on  stress-strain  curve  data. 

A  =  1.142  x  10'2  sec'1 

r  =  3.52 

Z.  =  100  ksi 
1 

This  second  group  of  constants  is  based  on  creep  test  data.  The  elastic 
modulus  at  1350°F  was  determined  to  be 

E  =  26.3  x  106  psi 

and  Poisson's  ratio  was  arbitrarily  chosen  as 

v  -  .3 
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Figure  4  shows  a  comparison  of  the  stress-strain  behavior  of  the 
Bodner  model  (using  the  given  IN-100  material  constants)  and  experimental 
stress-strain  data.  Each  of  the  curves  displays  the  response  of  the 
Bodner  model  under  the  loading  condition  of  a  given  constant  strain  rate. 

-3  -1 

The  experimental  data  were  generated  at  a  strain  rate  of  .83  x  10  sec 
and  compare  quite  well  with  the  Bodner  curve  for  the  same  strain  rate. 
Some  variation  occurs  in  the  region  of  initial  inelastic  behavior  but 
this  difference  is  small. 

Figures  5  and  6  compare  creep  behavior  of  the  Bodner  model  with 
experimental  data  at  the  stress  levels  of  127.3  ksi  and  72  ksi,  respec¬ 
tively.  Note  the  initial  experimental  curve's  slope  or  strain  rate 
is  duplicated  by  the  Bodner  model  in  both  figures.  However,  the  strain 
magnitudes  differ  somewhat  due  to  the  initial  time  required  for  the 
Bodner  model  to  reach  steady  state  creep  at  these  stress  levels.  Also 
as  the  experimental  strain  rate  increases  with  time  the  Bodner  model 
cannot  closely  follow  since  its  formulation  restricts  it  to  only 
secondary  type  creep  in  this  situation  (i.e.,  constant  strain  rate  for 
constant  stress). 

2.  SOLUTION  PROCEDURE  FOR  ELASTIC-VISCOPLASTIC  STRUCTURES 

In  elastic-plastic  analysis,  it  is  necessary  to  apply  loads 
incremental ly  to  satisfy  the  appropriate  yield  condition  (e.g.,  von 
Mises)  and  flow  rule  (e.g.,  Prandtl-Reuss )  associated  with  incremental 
plasticity  (Reference  64).  Similarly,  with  elastic-viscoplastic  behavior 
an  incremental  procedure  is  required,  but  here  time  is  incremented 
directly  while  load,  strain,  stress,  etc.  are  incremented  indirectly 
through  the  time  integration  procedure.  The  algorithm  used  for  a 
typical  time-step  in  the  elastic-viscoplastic  residual  force  method 
(Reference  10)  is  summarized  as  follows: 

1.  Add  time  increment  dt’  to  the  preceding  time  ti_1  to  obtain 

the  current  time  t1 . 
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Figure  4.  Stress-Strain  Behavior  of  Bodner  Model  (IN-100  Constants) 
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Figure  5.  Creep  Behavior  of  Bodner  Model  at  127.3  ksi  (IN-100  Constants) 
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Figure  6.  Creep  Behavior  of  Bodner  Model  at  72.0  ksi  (IN-100  Constants) 


Pi  *  P  i  1 

2.  Compute  increments  of  plastic  strain,  {de .  ^1  =  f e ^ ^ >  dt  and 

Pi  Pi-1  Pi 

add  to  preceding  plastic  strain,  {e.^}  =  +  *  *n 

general  the  plastic  strain  rate,  {e.J,  is  a  function  of  the 

1  J 

material's  viscosity  and  the  given  stress  level  (see  the  Viscoplastic 
Material  Models  section). 

3.  Compute  the  plastic  load  vector  {Q}^  =  [B3  [D]  {e^j}1  dvol . 

4.  Compute  the  current  external  load  vector  {P}’  =  {Pl^dt1  +  { P> ^  . 

5.  Compute  the  nodal  displacements  { U> 1  from  Equation  11, 

{U}1  =  [K]"1  {  {P}1  +  {Q}1"1}. 


6.  Compute  the  current  total  strain  {e. -}1  from  the  strain  dis- 

"  J 

placement  relationship,  { j > 1  =  [BjdJ}1. 

7.  Compute  the  current  stress  {o}1  as  follows,  {o^}1  = 

[D]{ 


8.  Check  the  time-step  size  in  terms  of  prescribed  stress  and 
strain  change  tolerances  per  time-step  (see  the  following  section 
on  Variable  Time  Step  Integration  of  Viscoplastic  Flow  Laws).  If 
these  tolerances  are  not  exceeded  the  time-step  size  may  be 
increased  for  the  next  time-step  or  left  the  same  value.  But 
if  the  tolerances  are  exceeded,  the  time-step  size  is  reduced 
and  steps  1  through  8  are  repeated  for  this  same  time-step  in 
an  effort  to  satisfy  the  stress  and  strain  change  tolerances. 


9.  Repeat  steps  1  through  8  until  the  desired  simulation  time 
is  reached. 
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3.  VARIABLE  TIME-STEP  INTEGRATION  OF  VISCOPLASTIC  FLOW  LAWS 


One  of  the  final  objectives  in  the  development  of  VISCO  was  to  be 
able  to  accurately  model  nonlinear  material  behavior  both  in  the  high 
strain  rate  load  up  stage  and  during  the  low  strain  rate  constant  load 
creep  stage.  For  stable  accurate  results,  the  time-step  size  must  be 
orders  of  magnitude  less  during  the  load  up  stage  compared  to  the  creep 
stage.  Therefore  to  be  computationally  feasible  in  a  large  problem 
(many  degrees  of  freedom)  some  method  is  necessary  to  determine  the 
maximum  time-step  during  each  stage  of  t.he  analysis  while  maintaining 
reasonable  accuracy.  Cormeau  (Reference  65)  investigated  the  numerical 
stability  of  simple  time  marching  schemes  used  in  elastic-viscoplastic 
analysis.  The  Malvern  and  Norton  models  were  studied  among  others. 

For  the  Malvern  model  the  following  maximum  time-step  size  was 
determined 


dtM  1 


4  ( 1 +v )  o 


(38) 


where  v  is  Poisson's  ratio  and  all  other  parameters  are  as  defined 
earlier.  For  the  Norton  model  the  maximum  time-step  was 


dtN  — 


4  (1+v) 

3YCEBoeB~1 


(39) 


In  general  dt^  is  several  orders  of  magnitude  larger  than  dt^.  If  the 

Malvern  model  or  the  Norton  model  are  used  separately  Cormeau  has  shown 
Equations  38  and  39  to  work  well.  However,  if  the  Norton  and  Malvern 
models  are  coupled  together  for  a  more  complete  flow  law  a  method  is 
necessary  to  provide  the  maximum  time-step  during  transition  from  load 
up  (Malvern  dominated  phase)  to  creep  (Norton  dominated  phase).  In 
addition,  this  time-step  maximizing  scheme  is  all  the  more  required 
when  numerically  integrating  the  Bodner  equations  since  Cormeau's 
analysis  does  not  directly  apply  to  the  Bodner  viscoplastic  flow  law. 
Consequently,  the  following  logic  which  is  also  similar  to  the  MARC 
program  (Reference  57),  was  incorporated  into  VISCO. 


3  c 
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This  time-step  maximizing  logic  continually  tries  to  increase  the 
time-step  size  subject  to  a  stress  and  strain  constraint.  These  con¬ 
straints  are  the  allowable  amounts  of  change  in  stress,  ato-| ,  and 

strain,  etoi ,  during  a  given  time-step  and  their  values  will  be 

discussed  later  in  this  section.  The  computer  algorithm  is  based 

around  parameters  P  and  P  defined  as  follows 

a  e 


P 

a 


P 

e 


— i  ■  ■ 

etotal  Etol 


(40) 


(41) 


where  e 


total 


V 


2  2  2 

e..  +  e_  +  . 50  y  and  the  superscript  i  refers  to 

xy 


'x  ~z 

the  time-step.  Note  that  if  the  change  in  effective  stress  between 
time-steps  i  and  i-1  just  satisfies  the  stress  constraint  or  stress 
tolerance,  ,  then  Equation  40  will  give  a  value  of  unity  for  P^. 

Similarly,  if  the  effective  plastic  strain  increment  for  time-step  i, 

P  i 

(deg)  ,  just  satisfies  the  strain  constraint  or  strain  tolerance,  otol 
then  P^  will  equal  unity  from  Equation  41. 


The  parameters  P  and  P  are  calculated  for  every  element  and  P 

c  e 

is  set  equal  to  the  largest  one.  One  method  of  changing  the  time-step 
size  based  on  P  is 

dt1  =  dt1-1/P  (42) 

Note  that  if  P  is  unity  no  change  in  the  time-step  size,  dt,  occurs. 
However,  if  P  is  less  than  one,  dt1  is  greater  than  dti_^  and  if  P  is 
greater  than  one  dt1  is  less  than  dti_V  In  the  case  of  P  being  greater 
than  one,  the  amount  of  change  in  stress  or  strain  has  exceeded  its 
respective  tolerance  and  recalculations  for  that  time-step  are  necessary 
using  the  reduced  time-step  size  from  Equation  42. 


AFWAL-TR- 80-41 40 


To  avoid  several  successive  recalculation  steps  that  develop  when 
P  is  greater  than  one  and  Equation  42  is  used,  the  following  substitute 
for  Equation  42  was  employed. 


dt1 

=  .8  dt1" Vp 

If 

P  >  1 

dt1' 

-  dt1'"1 

If 

.8  <  P  <  1 

dt1 

=  1.25  dt1"1 

If 

.65  <  P  <  .8 

dt1 

=  1.5  dt1"1 

If 

P  <  .65 

Note  that  Equation  43  reduces  the  time-step  size  more  than  Equation  42 
if  P  is  greater  than  one  but  if  P  is  less  than  one  Equation  43  increases 
the  time-step  size  slower  than  does  Equation  42.  Both  of  these  differ¬ 
ences  between  Equations  43  and  42  tend  to  reduce  the  number  of  calcula¬ 
tion  steps  required. 

Determination  of  the  values  for  the  stress  and  strain  tolerances, 
ot()1  and  etQ|  respectively,  was  accomplished  by  employing  VISCO  to 

analyze  a  plate  with  a  V-notch.  This  particular  geometry  was  chosen 
since  it  has  a  high  stress  concentration  around  the  notch  and  is  there¬ 
fore  somewhat  similar  to  a  plate  with  a  crack  which  is  the  geometry  of 
ultimate  interest  in  this  research  effort.  But  the  V-Notch  geometry 
can  be  modeled  with  far  less  elements  than  a  cracked  plate  requires  and 
still  compare  with  other  V-notch  solutions  in  the  literature  {Reference 
40).  Figure  7  shows  the  finite  element  mesh  employed  for  the  V-notch 
plate.  Only  or,e  qyarter  of  the  plate  is  modeled  due  to  symmetry. 

Element  sizes  were  made  smallest  near  the  V-notch  in  order  to  capture 
the  stress  concentration  there.  A  total  of  182  constant  strain  triangu¬ 
lar  elements  were  employed  which  is  somewhat  less  than  Yamada  et.  al . 
(Reference  40)  used.  However,  good  agreement  with  Yamada  was  achieved 
and  this  will  be  demonstrated  later  in  the  Validation  Examples  section. 
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Figure  7.  V-Notched  Plate  "4odel 
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a.  Malvern  Model 

—  p 

An  elastic-perfectly  VISCO  plastic  (i.e.  a(ee)  =  constant)  plane 

stress  analysis  was  carried  out  using  the  Malvern  model.  The  ratio  of 
elastic  modulus  E  to  the  yield  stress  was  288  and  the  ratio  of  applied 
remote  stress  to  yield  stress  was  0.593.  Poisson's  ratio  was  0.2.  The 
load  was  applied  in  nondimensional  time  of  v  t  =  10"^  where  y  is  the 

D  P 

fluidity  constant  in  the  Malvern  model  (Equation  14).  Results  were  taken 
after  all  stress  and  strain  rates  were  zero.  This  may  be  defined  as  the 
steady  state  condition  which  was  observed  to  occur  after  y  t  =  0.4. 

Several  analyses  of  the  V-notch  were  performed  using  different 
values  for  the  stress  and  strain  tolerances.  Examination  of  the  results 
found  that  between  displacements  near  the  notch,  plastic  strain  near 
the  notch,  and  total  plastic  strain  energy  in  the  model,  the  last  one 
was  most  sensitive  to  stress/strain  tolerance  variations.  Table  1 
displays  the  percent  variation  in  plastic  strain  energy  as  stress/strain 
tolerances  are  varied.  The  percent  variation  is  relative  to  the  plastic 
strain  energy  calculated  when  otol  =  .01  and  eto1  =  .01.  It  was  noted 

that  the  stress/strain  tolerances  that  kept  the  percent  variation  in 
plastic  strain  energy  around  one  percent  also  kept  the  time-step  size, 
during  most  of  the  computing,  under  Cormeau's  critical  value,  dtm> 
in  Equation  38.  The  amount  of  computer  time  increases  rapidly  as  the 
tolerances  are  reduced  to  a  value  of  .01.  A  good  compromise  between 
computer  time  required  and  apparent  accuracy  from  Table  1  is  a  stress 
tolerance  of  .03.  Note  that  for  this  stress  tolerance,  the  strain 
tolerance  can  be  relaxed  all  the  way  to  .20  with  little  change  in 
plastic  strain  energy.  These  results  are  in  agreement  with  the 
recommendations  for  stress/strain  tolerances  in  the  MARC  program 
(Reference  57). 

b.  Bodner  Model 

A  similar  stress/strain  tolerance  investigation  was  performed  with 
the  Bodner  material  model  in  plane  stress.  The  V-notch  model  in  Figure 


TABLE  1 
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11  was  again  employed.  The  material  constants  used  for  the  Bodner 
equations  were  those  matched  to  IN-100  and  given  previously.  The  ratio 
of  applied  remote  stress  to  yield  stress  was  again  .593  which  Is  based 
on  140  ksi  for  the  Initial  yield  stress.  Maximum  load  was  reached  in 
10  seconds  and  plastic  strain  energy  results  were  then  taken  similar 
to  the  Malvern  model  investigation.  Again  several  analyses  of  the 
V-notch  were  performed  using  different  combinations  of  stress  and 
strain  tolerance  values.  The  results  for  the  Bodner  model  are  given 
in  Table  2.  Comparing  Tables  1  and  2  shows  that  the  Bodner  model 
develops  less  plastic  strain  energy  variation  for  the  given  stress/ 
strain  tolerances  than  does  the  Malvern  model.  This  comparison  indicates 
that  higher  stress/strain  tolerances  could  then  be  used  for  the  Bodner 
model.  However,  if  the  analyses  are  continued  into  the  creep  phase, 
it  was  observed  that  the  stress  values  tend  to  oscillate  somewhat  when 
they  should  be  monotonically  relaxing  near  the  notch  tip.  This  oscilla¬ 
tion  was  fairly  well  damped  out  when  a  stress  tolerance  of  0.03  was  used. 
Once  this  small  value  for  the  stress  tolerance  is  chosen  the  strain 
tolerance  has  little  effect  as  long  as  it  is  greater  than  or  equal  to 
the  stress  tolerance. 

The  average  computation  time  for  the  Bodner  model  was  about  twice 
that  required  for  the  Malvern  model.  The  average  central  processor 
time  required  on  the  CDC-6600  was  150  seconds  for  this  200  degree  of 
freedom  problem  with  the  Bodner  model. 

4.  VALIDATION  EXAMPLES 

The  following  five  examples  were  performed  to  demonstrate  the 
validity  of  the  VISC0  computer  program.  The  first  four  examples  employ 
the  Malvern  model  to  compare  with  published  time  independent  elastic- 
plastic  solutions.  The  Bodner  model  is  then  tested  in  the  fifth  example 
by  comparing  its  behavior  with  the  results  from  coupling  the  Malvern 
and  Norton  flow  laws  together. 
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a.  Infinite  Sheet  with  Pressurized  Hole 


The  first  example  employs  the  Malvern  model  to  analyze  the  infinite 
sheet  with  a  pressurized  hole.  Comparisons  are  made  with  the  exact 
elastic-plastic  strain  history  independent  deformation  solution  by  Hsu 
and  Forman  (Reference  70).  The  material  properties  for  this  plane  stress 
solution  were  E  =  10.5  x  10^  psi,  v  =  0.5,  and  a  Ramberg-Osgood  stress- 
strain  equation 

e  =  for  lol  >  (44) 


where  e  is  the  total  strain  (i.e.,  elastic  plus  plastic),  o'  is  the 

initial  yield  stress  (o  =  55000  psi)  and  n  determines  the  shape  of  the 

curve  (n  =  9).  Equation  44  must  be  inverted  so  that  yield  stress  o  is  a 
function  of  plastic  strain  ep  in  order  to  be  used  in  the  Malvern  model. 
The  first  step  in  this  in  ersion  process  was  to  split  the  total  strain 
e  into  its  elastic,  eE,  and  plastic,  ep,  parts  and  then  set  the  elastic 
strain,  eE,  equal  to  o/E.  Equation  44  may  then  be  rewritten  as 


e 


(45) 


Solving  for  o11  results  in 

°n  ~  °o  (°  +  EeP)  (46) 

n+^ 

and  taking  the  n  root  of  Equation  46 

n-1  1 

a  =  <7q  n  (a  +  EsP)n  (47) 

The  stress,  a,  on  the  left-hand  side  of  Equation  47  will  now  be 

redefined  as  the  strain  hardening  yield  stress  o  so  that,  approximately 

n-1  D  i  _ 

a  =  oc  n  (o  +  Eek)  “  for  i a |  >  aQ  (48) 

Note  when  a  equals  o0  and  eP  is  zero,  Equation  48  is  identically  satisfied. 
Also  when  a  is  greater  than  a0,  the  plastic  strain  will  be  greater  than 
zero  and  the  product  E  e  P  will  have  the  major  effect  on  7,  which  is 
desired. 
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The  finite  element  mesh  for  modeling  the  infinite  plate  was  similar 
to  Figure  9b,  however,  132  triangular  elements  were  used.  The  outer 
radius  to  inner  radius  ratio  of  the  element  mesh  was  15  which  was 
assumed  to  approximate  an  infinite  radius  plate. 

Figure  8  shows  the  radial  and  tangential  stress  profiles  for  three 
internal  pressure  ratios,  P/o0.  Hsu  and  Forman  indicated  that,  in 
consideration  of  Budiansky's  criterion  for  the  acceptability  of 
deformation  theory,  their  infinite  plate  solution  should  not  disagree 
greatly  with  an  incremental  elastic-plastic  solution.  One  may  observe 
from  Figure  8  a  close  approximation  between  the  techniques  thus  lend¬ 
ing  validity  to  the  incremental  plane  stress  approach  incorporated 
within  VISCO. 

b.  Thick  Cylinder  with  Internal  Pressure 

The  second  example  employs  the  Malvern  model  to  analyze  a  thick 
cylinder  with  internal  pressure.  In  this  case  comparison  is  made  with 
a  non-finite  element  solution  by  Hodge  and  White  (Reference  71)  who 
again  used  deformation  theory.  This  was  an  elastic-perfectly  plastic 
plane  strain  analysis.  The  ratio  of  the  elastic  modulus  to  yield 
stress  was  E/a  =  190.9  with  v  =  0.33.  The  ratio  of  outer  radius  to 
inner  radius  for  the  cylinder  is  two.  Figure  9  shows  the  finite  ele¬ 
ment  mesh  used  to  model  a  symmetric  section  of  the  thick  cylinder. 

The  mesh  incorporates  40  triangular  elements.  Figure  9  shows  the 
tangential  stress  profile  for  a  pressure  ratio,  P/<7,  of  0.76.  Both 
the  elastic  and  the  elastic-plastic  profiles  are  given.  There  is 
good  agreement  with  Hodge  and  White  and  thus  validity  is  again  given 
to  the  incremental  plane  strain  portion  of  the  VISCO  program. 

c.  V-Notched  Plate  in  Tension 

The  third  example  employs  the  Malvern  model  to  analyze  a  V-notched 
plate  in  tension.  Comparison  is  made  with  another  finite  element  solu¬ 
tion  by  Yamada,  et.  al.  (Reference  40).  A  time  independent  tangent 
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Figure  9a.  Thick  Cylinder  with  Internal  Pressure,  Tangential  Stress 
Profile 
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Figure  9b.  Thick  Cylinder  with  Internal  Pressure,  Finite  Element  Mesh 
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modulus  approach  was  used  by  Yamada  et.  al.  for  their  elastic-plastic 
analysis.  This  was  an  el astic-perfectly  plastic  plane  stress  solution 
where  the  ratio  of  the  elastic  modulus  to  yield  stress  E/o  =  666.7  and 
v  =  0.3.  The  finite  element  mesh  used  herein  was  the  same  as  given 
previously  in  Figure  7  where  only  one  quarter  of  the  plate  is  modeled 
due  to  symmetry.  Element  sizes  were  made  smallest  near  the  V-notch 
in  order  to  capture  the  stress  concentration  there.  A  total  of  182 
constant  strain  triangular  elements  were  employed  which  is  somewhat 
less  than  used  by  Yamada  et.  al . 

Figure  10  shows  two  y-component  of  stress  profiles  for  the  minimum 
section  of  the  V-notched  plate.  Note  for  the  applied  stress  of  of  a  = 
0.185,  which  was  Yamada's  initial  yield  load,  the  present  results  from 
VISCO  agree  very  well  except  for  the  elements  immediately  at  the  notch. 
The  high  elastic  stress  concentration  near  the  notch  was  not  modeled 
as  well  by  the  present  VISCO  analysis  due  to  larger  element  sizes  being 
used  in  the  notch  vicinity.  However,  with  yielding  at  the  applied 
stress  of  afo  =  0.584,  comparison  with  Yamada's  results  at  the  notch 
are  much  better.  Note  that  plastic  action  diminishes  stress  gradients 
and  thus  fewer  elements  are  needed  for  an  elastic  solution  can  be  used 
to  produce  good  stress  profiles  after  yielding  occurs.  However,  it 
should  be  kept  in  mind  that  plastic  action  does  not  diminish  the  strain 
gradient  like  the  stress  gradient  and  therefore  strain  profiles,  even 
after  yielding,  will  be  quite  sensitive  to  element  sizing. 

Figure  11  displays  the  finite  element  mesh  for  the  V-notched  plate 
with  those  elements  left  out  that  have  exceeded  the  yield  stress  for  the 
applied  stress  of  a/a  =  0.584.  The  absence  of  these  elements  thus 
describes  the  plastic  zone  size  and  compares  well  with  that  of  Yamada. 
Numbers  within  the  elements  of  Figure  11  reflect  each  elements  effective 
stress  as  a  percentage  of  yield  stress. 
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Figure  10.  Stress  Profile  for  V-Notched  Plate  in  Tension 
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Figure  11.  Plastic  Zone  in  V-Notched  Plate 
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d.  Cracked  Three  Point  Bend  Specimen 

The  fourth  example  employs  the  Malvern  model  to  analyze  a  cracked 
three-point  bend  specimen  as  shown  in  Figure  12.  This  will  be  the  same 
problem  as  was  solved  with  ten  different  computer  codes  for  an  ASTM 
analytical  round  robin  and  documented  by  Wilson  and  Osias  (Reference 
72).  This  was  a  plane  strain  crack  problem  where  E=31 . 59  x  106  psi, 
v=  0.30,  and  the  following  Ramberg-Osgood  equation  was  used  for  the 
stress-strain  curve 


where  e  is  the  total  strain  (i.e.,  elastic  and  plastic),  n  is  10,  and 
Bg  is  120  x  10^  psi.  To  make  Equation  49  compatible  with  the  Malvern 

model,  stress  must  be  written  as  a  function  of  plastic  strain.  There¬ 
fore,  after  subtracting  the  elastic  strain  (also  equal  to  a/E)  from  both 
sides  of  Equation  49,  solving  for  the  stress  o,  and  redefining  o  as  a 

1 

a  »  Bg  (eP)  n  (50) 

Due  to  symmetry  only  one  half  of  the  three-point  bend  specimen  was 
modeled  with  the  finite  element  mesh  in  Figure  13.  This  particular 
pattern  for  the  element  mesh  was  used  by  Ohtani  and  Nakamura  (Reference 
61)  for  a  center  cracked  plate.  Note  how  the  element  sizes  are  reduced 
as  indicated  in  Figure  13  by  arrows  to  the  first  and  second  reduction. 
This  pattern  provides  for  an  unlimited  number  of  element  size  reductions 
while  also  maintaining  good  element  aspect  ratios  (e.g.,  from  1  to  0.5) 
and  ensuring  no  two  neighboring  elements  differ  in  size  by  more  than  a 
factor  of  2.  Accordingly,  each  reduction  cuts  the  preceding  element 
size  in  half.  Eight  element  size  reductions,  incorporated  in  the  Figure 
13  mesh,  stepwise  reduced  the  0.20  inch  elements  at  the  upper  boundary 
to  a  7.8125  x  10”^  inch  element  at  the  crack  tip.  This  crack  tip  element 
size  was  slightly  smaller  than  the  smallest  used  in  Reference  72.  The 
total  number  of  elements  in  Figure  13  was  584  with  388  nodes.  Figures 
14,  15,  and  16  show  the  results  of  V I  SCO  compared  to  ten  times  indepen¬ 
dent  elastic-plastic  solutions  documented  by  Wilson  and  Osias.  Eight  of 
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Figure  12.  Three-Point  Bend  Specimen  with  Crack 
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Figure  14.  Crack  Opening  Displacement  Profile  @  UA  =  .050  In.  for 
Three-Point  Bend  Model 


AFWAL-TR-80-41 40 


these  solutions  fall  within  lines  II  and  III  in  Figures  15  and  16.  The 
other  two  solutions  fall  near  line  I  in  these  two  figures.  In  all  cases 
the  present  results  from  VISCO  fall  within  the  ASTM  round  robin  range 
as  shown  in  these  figures.  Therefore  the  element  sizing  and  arrangement 
in  Figure  13  as  used  herein  with  VISCO,  appears  to  give  good  results  for 
a  complex  problem  that  includes  a  crack  and  load  that  induces  bending. 


e.  Center  Cracked  Plate  with  Bodner  Model 


The  fifth  example  employs  the  Bodner  model  to  analyze  a  center 
cracked  plate.  Comparison  will  be  made  with  results  from  a  similar 
analysis  using  the  Malvern  model  coupled  with  Norton's  Law  as  follows 


•  P 

eij 


Cy. 


fer  - ') 


+  Yc  {ae}  ]  f  ^  if  °e  "  a{tl] 


Cv, 


(°je]  f  ^  if  °e  i 


(51 


This  Mai vern-Norton  combination  is  a  superposition  approach  suggested 
by  Zienkiewicz  (Reference  10)  to  model  in  a  unified  sense  both  initial 
load  up  viscoplasticity  and  creep  under  sustained  load.  Also,  to  assess 
the  contribution  of  pure  secondary  creep  in  the  Bodner  model,  comparison 
will  be  made  to  an  analysis  using  only  Norton’s  Law  for  the  viscoplastic 
material  model. 


The  material  properties  will  be  those  matched  to  IN-100  at  1350°F. 
The  constants  for  the  Bodner  model  will  be  those  previously  given  for 
IN-100.  The  yield  stress  as  a  function  of  plastic  strain,  oU^), 
needed  for  the  Malvern  model  will  be  a  multilinear  fit  to  the  experi¬ 
mental  stress-strain  curve  in  Figure  4  and  given  as  follows 


/ 
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The  fluidity  constant  yp  for  the  Malvern  model  will  be  given  a  value 

of  .58  sec”^ .  This  value  for  yp  is  sufficiently  high  such  that  the 

given  o’  function  is  followed  very  closely  (e.g.,  recall  Figure  3  for 
Yp  =  «).  The  constants  for  Norton's  Law  will  be  yc  =  3.7394  x  10"®^ 

(psi )~l°'64/sec  and  g  =  10.64.  These  values  were  determined  by  matching 
the  initial  IN-100  creep  behavior  in  Figures  5  and  6  as  discussed  in  the 
Norton's  Law  for  secondary  creep  section. 

Due  to  symmetry  only  one  fourth  of  the  center  cracked  plate  is 

modeled  by  the  finite  element  mesh  given  in  Figure  17.  The  plate's 

height  is  2.8  inch,  width  is  1  inch,  and  thickness  is  0.3  inch.  The 

crack  length,  2a,  is  .2734  inch  or  a/W  equal  to  .1367.  The  numbers 

inside  the  elements  indicate  a  total  of  355  elements  were  employed. 

The  total  number  of  nodes  is  211.  The  triangular  elements  around  the 

-4 

crack  tip  have  a  height  and  base  dimension  of  7.8125  x  10  inch. 

Further  discussion  of  this  element  mesh  will  be  presented  in  the 
Applications  section  since  this  particular  mesh  was  also  employed  there 
to  simulate  the  experimental  program.  A  maximum  stress  of  36320  psi 
will  be  applied  to  the  upper  boundary  of  the  center  cracked  plate  in 
five  seconds  and  then  held  constant  for  an  elapsed  time  of  1000  seconds. 
This  applied  stress  level  corresponds  to  a  load  level  also  used  in  the 
experimental  program  to  be  simulated  later. 

A  comparison  of  the  behavior  of  these  three  material  models  is 
given  in  Figures  18,  19,  and  20.  The  effective  stress  and  plastic 
strain  in  Figures  18  and  19  are  from  the  element  at  the  crack  tip 
which  had  the  highest  elastic  stress  concentration  factor. 

The  stress-strain  behavior  in  Figure  18  occurred  over  a  time  period 
of  1000  seconds.  The  values  to  the  left  of  6%  strain  occurred  in  approxi¬ 
mately  five  seconds  (the  load  up  period)  whereas  to  the  right  of  6% 
strain,  stress  relaxation  and  redistribution  is  taking  place  over 
approximately  1000  seconds  of  sustained  load  creep  behavior.  Note  how 
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Figure  17.  Center  Cracked  Plate  Finite  Element  Model,  a/W  =  0.1367 


0.7  in. 


EFFECTIVE  STRESS  (ksi) 


AFWAL-TR-80-4140 


18)  f 


I 


160 


L  A- 
& 


m 


w 


o  o 


— Norton’s  Model 

\ 

G 

O 

0  N 

A 

A 

A  0A 

n 


o  Malvern-Norton  Model 
a  Bodner  Model 
□  Norton  Model 


2  H  6  8  10 

EFFECTIVE  PLASTIC  STRAIN  (%) 


12 


Figure  18.  Stress  vs.  Strain  at  the  Crack  Tip  for  Bodner,  Malvern-Norton, 
and  Norton  Material  Models 
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Figure  20.  Crack  Mouth  Displacements  in  Center  Cracked  Plate  for  Bodner, 
Mai vern-Norton,  and  Norton  Material  Models 
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the  Norton  model  behaves  during  the  load  up  period  (i.e.,  strains  less 
than  6%).  Due  to  its  slow  creep  response,  the  effective  stress  from  the 
Norton  model  approaches  the  values  that  would  occur  during  load  up  in  an 
elastic  analysis.  However,  the  Norton  model  relaxes  the  stress  down  to 
values  very  similar  to  the  other  material  model  stress  values.  In 
general  the  Bodner  model  behaves  very  similar  to  the  Mai vern-Norton 
combination  both  during  load  up  and  the  sustained  load  creep  periods. 

Figure  19  shows  the  time  dependent  behavior  of  the  effective  stress 
and  plastic  strain  from  the  same  crack  tip  element.  After  approximately 
200  seconds,  all  three  material  models  develop  nearly  identical  effective 
stress  values.  Plastic  strain  behavior  with  time  is  also  very  similar 
for  the  three  models  except  for  under  100  seconds  of  time.  The  differ¬ 
ence  in  percent  plastic  strain  between  the  three  curves  developed  pri¬ 
marily  from  different  plastic  strain  rates  during  the  load  up  phase  and 
remains  fairly  constant  for  time  greater  than  200  seconds. 

Figure  20  shows  how  the  crack  mouth  displacement  increases  with  time 
after  the  maximum  load  is  achieved.  The  location  of  the  crack  mouth 
displacement  is  indicated  in  Figure  20  to  be  0.050  inches  from  the 
vertical  centerline.  Again  the  three  curves  are  very  similar  after  200 
seconds  and  their  separation  is  due  primarily  to  dissimilar  behavior  for 
time  under  200  seconds.  The  Norton  model  displays  more  displacement 
after  maximum  load  than  the  other  models  since  it  is  effectively  making 
up  for  its  slower  plastic  strain  rates  and  associated  displacements 
during  the  load  up  phase.  This  apparently  is  also  true  when  comparing 
the  Mai vern-Norton  to  the  Bodner  model  displacement  curve,  however, 
here  the  difference  is  much  less. 

Therefore,  based  on  these  comparisons  with  somewhat  similar  material 
models  used  to  analyze  the  center  cracked  plate,  the  Bodner  model  within 
VISC0  is  considered  to  be  working  well.  For  further  discussion  and 
comparison  of  these  material  models  see  Reference  77. 
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5.  SUMMARY  OF  FINDINGS 

The  following  findings  were  determined  from  the  preceding  validation 
examples  and  the  general  development  of  the  VISCO  program. 

1.  The  stress  tolerance,  ,  which  controls  the  variable  time 

step  size  in  VISCO's  numerical  time  integration  algorithm, 
provides  good  results  for  reasonable  computer  time  requirements 
when  set  to  the  value  of  0.03.  The  strain  tolerance,  etQi ,  has 

little  effect  as  long  as  its  greater  than  or  equal  to  the 
stress  tolerance.  This  finding  pertains  to  the  Bodner,  Malvern, 
and  Norton  material  models. 

2.  The  VISCO  program,  while  employing  the  Malvern  material  model 
agrees  well  with  so-called  exact  elastic-plastic  deformation 
solutions  both  in  plane  stress  and  plane  strain  conditions. 
Agreement  is  also  good  with  time  independent  elastic-plastic 
finite  element  solutions  in  a  V-notched  plate  and  a  cracked 
three-point  bend  specimen  whose  results  came  from  an  ASTM 
analytical  round  robin. 

3.  The  finite  element  mesh  pattern  in  Figure  13  works  well  for 
modeling  cracked  plates.  This  pattern  conveniently  provides 
for  an  unlimited  number  of  element  size  reductions  to  capture 
the  crack  tip  singularity  while  also  maintaining  good  element 
aspect  ratios  and  minimizing  the  total  number  of  elements 
required. 

-4 

4.  A  crack  tip  element  size  to  specimen  width  ratio  of  7.8125  x  10 
in  a  cracked  three-point  bend  specimen  provided  good  agreement 
with  an  ASTM  analytical  round  robin. 

5.  The  Bodner  material  model  has  been  shown  to  behave  similarly  to 
the  Malvern-Norton  model  for  both  load  up  and  sustained  load 
creep  stages.  Also,  for  time  greater  than  200  seconds  after 
load  is  applied,  the  pure  secondary  creep  Norton  law  behaves 
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similar  to  the  Bodner  model.  In  general,  these  similarities  in 
material  model  behavior  should  be  true  for  most  any  metal,  how¬ 
ever,  the  indicated  200  second  time  delay  between  the  Norton 
and  Bodner  models  pertains  specifically  to  the  IN-100  alloy  at 
1 350°F. 

Therefore,  these  findings  support  the  validation  of  the  VISCO 
computer  program  and  provides  some  of  the  required  details  for  applying 
the  VISCO  program  with  the  Bodner  model  to  creep  crack  growth 
simulation. 
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SECTION  IV 

HYBRID  EXPERIMENTAL-NUMERICAL  PROCEDURE  TO 
ANALYZE  CREEP  CRACK  GROWTH 

The  simultaneous  use  of  experimental  data  from  crack  growth  tests 
and  a  theoretical  model  of  the  experimental  cracked  specimen  has  been 
termed  the  Hybrid  Experimental-Numerical  procedure  by  Kobayashi  (Refer¬ 
ence  16).  One  example  of  this  procedure  would  be  to  grow  a  crack  at 
experimentally  determined  rates  through  the  theoretical  model  of  the 
experimental  specimen.  Then*  from  the  results,  one  could  seek  out 
potential  crack  growth  criteria  which  hopefully,  during  crack  growth, 
display  themselves  as  fairly  constant  parameters  independent  of  both 
crack  length  and  load  intensity  (e.g.,  stress  or  strain  at  the  crack 
tip,  crack-opening  displacement,  etc.). 

In  the  present  investigation,  the  Hybrid  Experimental-Numerical 
procedure  required  the  theoretical  model  to  track  experimental  displace¬ 
ment  rates  rather  than  crack  growth  rates.  The  theoretical  model 
consisted  of  the  VISCO  program  employing  the  Bodner  constitutive 
equations  and  a  finite  element  mesh  of  the  experimental  specimen.  The 
experimental  displacements  were  either  measured  near  the  crack  tip  or 
near  the  vertical  centerline  (crack  mouth)  of  the  center  cracked  plate 
test  specimens  as  shown  in  Figure  21.  The  displacements  were  measured 
continuously  with  time  by  an  optical  interferometric  displacement 
measurement  technique  developed  by  Sharpe  (Reference  74).  Figure  22 
shows  a  typical  experimental  displacement  versus  time  curve.  The 
present  numerical  procedure  required  the  crack  to  grow  through  the 
theoretical  model  such  that  the  displacements  due  to  elastic-plastic 
behavior  and  crack  growth  added  up  to  the  experimental  displacements 
as  time  progressed.  Therefore,  crack  length  versus  time  became  a 
product  of  the  present  analysis  rather  than  an  input. 

1 .  CRACK  LENGTH  VERSUS  TIME 

Ideally  for  this  investigation,  the  experimental  data  should  be  in 
the  form  of  crack  length  versus  time  rather  than  displacement  rates. 
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Figure  21.  Center  Cracked  Plate  Test  Specimen 
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Figure  22.  Displacement  Under  Constant  Load 
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Unfortunately,  it  was  nearly  impossible,  experimentally,  to  measure  crack 
length  as  a  function  of  time  with  any  degree  of  precision  or  reliability 
on  the  surface  since  the  total  amount  of  crack  growth  in  these  tests  was 
extremely  small  (e.g.,  100  microns)  and  since  the  creep  crack  would  grow 
internally  or  tunnel  without  associated  surface  crack  growth. 

Another  attempt  to  measure  crack  growth  rate  indirectly  employed 
the  elastic  compliance  method  first  demonstrated  by  Clarke  (Reference 
75).  This  method  utilizes  the  change  in  elastic  compliance  of  the 
specimen  with  time  and  then  with  the  aid  of  a  compliance  versus  crack 
length  relationship  based  on  linear  elastic  fracture  mechanics,  crack 
growth  with  time  may  be  determined.  Figure  23  shows  schematically  an 
experimental  creep  crack  growth  load  versus  time  history  designed  to 
provide  discrete  values  of  compliance  at  selected  time  intevals.  The 
load  is  reduced  approximately  20%  and  then  restored  as  shown  at  times 
t^  through  t3  to  provide  load  displacement  data  at  various  times  during 

the  test.  Figure  24  shows  a  typical  set  of  load  displacement  data  for 
increasing  times  t-j  to  t4  in  a  creep  crack  growth  test.  Compliance  is 

defined  as  displacement  divided  by  load  and  therefore  the  slope  of  each 
line  in  Figure  24  represents  the  compliance  at  each  particular  time. 

Note  that  compliance  is  shown  to  decrease  in  going  from  time  t-j  to  t^ 

and  then  increase  from  time  t^  to  t^.  Comparing  this  behavior  to  a 

typical  elastic  compliance  versus  crack  length  curve  given  in  Figure  25, 
the  mathematical  implication  is  that  the  crack  shortens  with  time. 
Although  this  compliance  decrease/increase  behavior  has  also  been 
observed  by  Donat  (Reference  76),  no  known  experimental  data  supports 
any  physical  shortening  or  healing  of  the  crack.  However,  the  important 
implication  of  this  compliance  behavior  is  that  a  one  to  one  relation 
between  crack  length  and  compliance  does  not  exist  during  creep  crack 
growth.  Therefore,  the  elastic  compliance  method  cannot  be  used 
directly  to  measure  crack  length  with  time  in  creep  crack  growth  tests. 
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Figure  24.  Experimental  Compliance  Changes  with  Time 
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Figure  25.  Elastic  Compliance  vs.  Crack  Length 
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2.  MATCHING  EXPERIMENTAL  DISPLACEMENT  DATA 

In  each  Hybrid  Experimental-Numerical  application  of  the  VISCO 
finite  element  model  to  simulate  creep  crack  growth  tests,  the  y-dis- 
placement  of  the  node  which  was  closest  to  the  point  where  the  experi¬ 
mental  measurement  was  made  was  monitored  with  time.  If  the  node's 
displacement  became  less  than  experimental  displacement  at  a  given 
time  a  crack  tip  node  would  begin  to  be  released  to  simulate  crack  growth 
through  the  model . 

3.  CRACK  TIP  NODE  RELEASE  METHODS 

Crack  tip  nodes  were  released  in  one  of  two  different  methods.  The 
first  method  releases  the  node  and  totally  unloads  it  in  five  seconds. 

The  second  method  unloads  the  current  crack  tip  node  linearly  with  time 
over  the  total  time  span  between  the  time  when  the  current  crack  tip 
node  is  begun  to  be  released  and  the  time  when  the  next  crack  tip  node 
will  be  released.  The  five-second  node  release  method  for  crack  growth 
must  be  used  when  matching  experimental  displacements  or  when  a  crack 
growth  criterion  is  used.  In  both  of  these  cases  when  certain  conditions 
are  satisfied,  the  crack  must  grow  so  a  node  is  released.  However, 
when  the  current  node  is  being  unloaded  it  is  not  known  when  the  next 
crack  tip  node  will  be  released  and  therefore  the  continuous  unload 
method  cannot  be  used.  The  five-second  unload  time  for  the  first  method 
was  based  on  the  size  of  the  crack  tip  elements  and  the  maximum  crack 
growth  rates  occurring  in  the  creep  crack  growth  test  data  (i.e., 
maximum  crack  growth  rate  equals  element  size  divided  by  five  seconds). 

If  a  crack  growth  rate  criterion  were  used  then  it  could  be  determined 
by  extrapolation  when  the  next  crack  tip  will  be  released  and  thus  the 
current  crack  tip  node  could  be  unloaded  in  a  continuous  fashion  by  the 
second  method.  The  second  node  release  method  can  also  be  employed 
when  all  node  release  times  are  specified  at  the  beginning  of  the 
computer  run  (e.g.,  release  times  based  on  results  from  a  prior  computer 
run  using  release  method  one  and  matching  experimental  data). 
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In  both  node  release  methods,  the  node  force  required  to  hold  the 
crack  tip  node  at  zero  displacement  would  be  calculated  from  the  stresses 
in  the  adjacent  elements  as  follows 

(f)  =  53  I  W  {oij}  d  vo1  (53) 

Adjacent  vol 
Elements 


which  is  consistent  with  the  formulation  of  the  stiffness  matrix  in 
Appendix  A.  The  crack  tip  node  restraint  force,  f  ,  is  then  the 


component  of  {f}  perpendicular  to  the  crackline.  The  boundary  condition 
on  the  node  is  converted  from  zero  displacement  to  a  force  equal  to  f  . 

This  force  f  is  then  removed  depending  upon  which  node  release  method 

is  chosen  (Figure  26). 


The  change  of  the  crack  tip  node's  boundary  condition  from  dis¬ 
placement  to  a  force  boundary  condition  is  handled  very  conveniently 
with  the  Gauss-Seidel  iterative  linear  equation  solver  as  discussed  in 
Appendix  B.  Whenever  a  node  is  fixed  in  a  certain  direction,  its 
equilibrium  equation  in  that  direction  is  skipped  over  during  the 
iterative  solution  procedure  and  when  the  node  is  released  its 
equilibrium  equation  is  included  within  the  iterative  procedure. 

No  lengthy  refactorization  of  the  stiffness  matrix  is  required  for 
these  boundary  condition  changes. 
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Figure  26a.  Crack  Tip  Node  Unload  Methods,  Finite  Element  Crack  Tip. 


Figure  26b.  Crack  Tip  Node  Unload  Methods,  Five  Seconds  Node  Unload  Method 


TIME 


Figure  26c.  Crack  Tip  Node  Unload  Methods,  Continuous  Node  Unload  Method 
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SECTION  V 

APPLICATION  OF  THE  HYBRID  EXPERIMENTAL 
NUMERICAL  PROCEDURE  TO  CREEP  CRACK  GROWTH 

The  first  objective  in  the  present  research  was  to  develop  a  method 
for  getting  crack  growth  behavior  solely  from  displacement  measurements 
made  on  a  cracked  specimen  under  constant  load  and  at  elevated  tempera¬ 
ture.  This  objective  is  extremely  important  since  small  but  significant 
amounts  of  crack  growth  can  not  otherwise  be  resolved  by  conventional 
experimental  crack  measuring  techniques.  The  second  objective  was  to 
seek  out  crack  growth  criteria  based  on  the  crack  growth  behavior  identi¬ 
fied  from  the  present  work  on  the  first  objective  and  by  examining  various 
parameters  around  the  crack  tip  in  the  theoretical  model  (e.g.,  stress, 
strain,  crack  opening  displacement,  etc.). 

This  section  presents  the  results  of  applying  the  hybrid  experi¬ 
mental  numerical  (HEN)  procedure  to  creep  crack  growth  in  IN-100  at 
1 350°F.  The  experimental  portion  of  the  procedure  consisted  of 
displacement  versus  load  (i.e.,  compliance)  and  displacement  versus  time 
test  data  reported  by  Sharpe  (Reference  15).  The  numerical  portion  of 
the  HEN  procedure  consisted  of  the  VISCO  finite  element  program 
employing  the  Bodner  material  model.  The  material  constants  for  the 
Bodner  model  were  those  given  in. Section  III. 

The  machining  specifications  for  the  specimens  used  in  the  experi¬ 
mental  program  are  shown  in  Figure  27.  Only  the  center  uniform  cross 
section  part  of  the  specimen  was  considered  in  the  VISCO  simulation  and 
due  to  symmetry  only  one  quadrant  of  this  section  was  represented  by 
the  finite  element  meshes  given  in  Figures  28,  29,  and  30.  Each  of 
these  meshes  represent  the  center  cracked  plate  test  specimen  with 
different  crack  lengths.  The  convergence  of  these  meshes  has  been 
verified  through  the  work  in  Section  III  and  further  discussion  of 
their  accuracy  will  be  provided  in  subsequent  paragraphs.  Figures 
28,  29,  and  30  have  half  crack  lengths,  a,  or  0.137  inches,  0.237  inches, 
0.312  inches,  respectively.  Due  to  the  variations  of  the  initial  crack 
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Figure  27.  Experimental  Creep  Crack  Growth  Specimen 
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a  =  .1367  in.  Cracktip 


W/2  =  0.5  in. 

Figure  28.  Center  Cracked  Plate  Finite  Element  Model,  a/W  =  0.1367 
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Figure  30.  Center  Cracked  Plate  Finite  Element  Model,  a/W  =  0.3117 
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length  through  the  thickness,  the  surface  measurement  cannot  be  used 
directly  and  thus  effective  initial  crack  lengths  must  be  determined. 
These  effective  initial  crack  lengths  were  determined  by  matching 
experimental  load-displacement  data  along  the  crack  such  that  the 
finite  element  model  displayed  elastically  the  same  compliance  as 
the  experimental  specimens  did.  For  efficiency  sake  it  has  been  found 
realistic  to  use  the  same  element  mesh  for  slightly  different  initial 
crack  lengths. 

The  VISCO  finite  element  program  has  both  plane  stress  and  plane 
strain  analysis  capability.  The  current  investigation  chose  plane 
stress  as  reported  in  Reference  78  where  theoretical  plane  stress  J 
integral  values  agreed  best  with  test  data  for  an  even  thicker  (  1  inch 
thick)  compact  tension  specimen.  A  theoretical  model  must  display 
realistic  compliance  behavior  in  order  to  calculate  J  values  that 
agree  with  test  data.  Likewise,  in  the  present  research  realistic 
compliance  behavior  is  a  necessity. 

Figures  28,  29,  and  30  display  the  same  general  pattern  of  elements 
which  worked  well  for  the  three-point  bend  specimen  in  Section  III.  The 
elements  at  the  crack  tip  have  a  height  and  width  of  7.8125  x  10"4  inch. 
This  size  crack  tip  element  in  combination  with  the  given  general  ele¬ 
ment  pattern,  provided  for  355,  278,  and  362  elements  respectively  in 
the  three  figures.  Figures  28  and  30  provided  20  uniform  elements  ahead 
of  the  crack  tip  for  subsequent  crack  growth  whereas  Figure  29  had 
eight  uniform  elements.  Figure  31  shows  the  expanded  element  layout 
around  the  crack  tip  from  Figure  29.  This  region  of  uniform  elements 
ahead  of  the  initial  crack  tip  avoids  unrealistic  changes  in  compliance, 
as  the  crack  grows  through  the  model,  that  can  develop  when  nonuniform 
element  sizes  are  used.  The  number  of  uniform  elements  incorporated 
ahead  of  the  crack  tip  is  a  compromise  between  anticipated  crack 
growth  and  computer  time  required. 
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The  elastic  compliance  of  the  specimens  modeled  by  these  three 
meshes  in  VISCO  was  compared  to  an  empirical  solution  by  Eftls  and 
Liebowitz  (Reference  79)  in  the  form  of 


-1 

cosh 


£s  ec(gi)] 


(54) 


This  comparison  is  shown  in  Figure  32,  with  the  symbols  a,v,  and  W  as 
defined  in  the  figure.  Also  shown  is  nondimensionalized  compliance 
versus  crack  length  for  locations  off  the  vertical  centerline  and 
relatively  near  the  crack  tip.  These  locations  correspond  to  the 
displacement  measurement  locations  used  in  the  experimental  program. 

It  can  be  seen  that  the  centerline  compliance  agrees  very  well  with 
published  results.  Also  note  that  if  the  off  centerline  results  are 
linearly  extrapolated  (i.e.,  dashed  line)  back  to  the  Effis  and  Liebowitz 
curve  they  intersect  at  a/W  values  which  correspond  to  their  distance 
behind  the  crack  tip.  These  results  all  support  the  validity  of  the 
finite  element  meshes  used  in  the  present  investigation. 


Figure  33  shows  elastic  crack  opening  displacement  profiles  using 
VISCO  and  the  mesh  in  Figure  28  compared  to  the  Westergaard  equation 
for  elastic  displacements  around  the  crack  tip.  The  following  form  of 
the  Westergaard  equation  is  restricted  to  plane  stress  displacements 
behind  the  crack  tip  and  on  the  crack  surface  (i.e.,  crack  opening 
displacement) 


v 


/T 

2ti 


(55) 


where  r  is  the  distance  behind  the  crack  tip  and  Kj  is  the  mode  I 

elastic  stress  intensity  factor  which  for  the  center  crack  plate  can 
be  written  as  (Reference  73) 


Kj  =  0/(ira)sec(*a.) 


(56) 


where  again  a  and  W  are  the  half  crack  length  and  specimen  width 
respectively.  Note  that  agreement  with  the  Westergaard  equation  is 
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Figure  32.  Comparison  of  Compliance  from  VISCO  with  Eftis  &  Liebowitz 
for  the  Center  Cracked  Plate 
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quite  good  here  for  the  results  from  VISCO  using  the  mesh  with  the 
shortest  crack  length  (Figure  28).  The  same  crack  tip  element  size 
was  used  in  all  three  crack  length  models.  Therefore,  agreement  with 
Equation  55  would  even  be  better  for  the  longer  crack  length  models 
since  the  ratio  of  element  size  to  crack  length  will  be  smaller.  It 
has  been  shown  through  elastic  finite  element  convergence  studies  that 
as  this  ratio  of  crack  tip  element  size  to  crack  length  decreases, 
accuracy  increases  (Reference  81). 

Another  consideration  supporting  this  size  crack  tip  element 
(i.e.,  7.8  x  10-4  in.)  is  that  IN-100  grains  are  approximately  the 
same  dimension  (Reference  13).  It  may  be  argued  that  the  finite 
element  method  is  a  continuum  analysis  tool  and  that  accuracy  should 
only  improve  as  elements  are  refined.  However,  realistically  a 
continuum  does  not  exist  at  and  below  the  grain  size,  especially 
around  the  crack  tip.  Thus  incorporation  of  elements  smaller  than 
the  grain  size  might  be  unrealistic.  Furthermore,  it  is  postulated 
that  the  physical  blunting  of  the  crack  tip  that  can  be  associated 
with  noncontinuous  grain  structured  material  might  be  more  effectively 
considered  by  the  finite  element  model  used  herein  since  the  crack 
tip  elements  have  a  dimension  on  the  order  of  a  grain  size. 

1.  CRACK  GROWTH  PREDICTIONS 

A  summary  of  the  creep  crack  growth  test  data  reported  by  Sharpe 
(Reference  15)  is  given  in  Tables  3  and  4.  In  general  the  experimental 
program  objectives  were  to  generate  creep  crack  growth  data  in  the 
center  cracked  plate  specimen  at  1350°F  for  several  different  loads 
and  cracked  lengths.  The  loads  were  specified  in  terms  of  a  range  of 
stress  intensity  factor  values  from  approximately  15.0  to  35.0  ksi  /Tn 
for  each  respective  initial  crack  length.  With  these  objectives  in 
mind  and  the  limited  number  of  test  specimens,  only  one  test  was  done 
at  each  of  the  test  conditions.  In  the  process  of  developing  the 
experimental  procedure,  several  tests  did  not  result  in  good  data  and 
consequently  were  not  used  in  the  presc.t  investigation  as  Implied 
by  the  discontinuous  test  numbers  in  Table  3. 
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The  initial  surface  crack  lengths,  a$,  were  measured  on  the  surface 
of  the  test  specimens  after  fatigue  precracking  and  prior  to  the  load 
application  for  the  creep  crack  growth  tests.  Displacements  across  the 
crack  were  meisured  at  the  indent  locations,  qualitatively  referred  to 
as  mouth  or  tip  locations  described  in  Section  IV.  In  each  of  these 
creep  crack  growth  tests,  the  load  was  applied  in  five  seconds  and  held 
constant  for  the  given  test  duration  time  (excluding  small  unload/reload 
cycles  for  compliance  referred  to  in  Section  IV). 

Figure  34  is  a  photograph  of  the  fracture  surface  of  specimen  number 
2.  This  photograph  was  taken  after  breaking  open  the  center  cracked  plate 
specimen  following  creep  crack  growth  tests  5  and  6.  The  two  dark  bands 
indicated  with  arrows  are  the  creep  crack  extensions  from  tests  5  and  6 
whereas  the  remaining  fracture  surface  is  either  pre-test  or  post- test 
fatigue  crack  growth. 

Post- test  crack  length  measurements  were  made  on  the  fracture  sur¬ 
face  as  described  in  Figure  35.  The  average  of  the  four  crack  length 
measurements  was  defined  as  the  experimental  crack  length  a  .  This 
averaging  was  done  to  smooth  out  crack  length  differences  due  to 
asymmetric  crack  growth  and  variations  through  the  thickness.  Table  4 
gives  a  tabulation  of  the  initial  crack  length,  ag  and  am.  In  addition, 
a  crack  length  determined  by  compliance,  ac>  and  the  initial  crack 
length  used  in  the  finite  element  model,  aQ,  is  given.  The  a£  crack 
length  was  determined  as  discussed  earlier  by  varying  the  crack  length 
in  the  finite  element  model  until  the  model's  compliance  matched 
experimental  compliance  data. 

For  convenience  a  modified  form  of  the  above  compliance  technique 
for  crack  length  determination  was  also  employed  which  reduced  the  num¬ 
ber  of  element  meshes  for  different  crack  lengths  required.  This 
modification  to  the  compliance  method  made  use  of  the  experimental 
compliance,  CE,  and  the  finite  element  model  compliance,  Cp^,  for  a 

crack  length  near  the  test  value.  To  determine  the  test  specimen's 
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crack  length,  the  difference  between  Cp  and  Cpp  was  divided  by  the  rate 

of  change  of  compliance  with  crack  length  Ac/Aa  (i.e.,  slope  of  curve 
In  Figure  32  as  follows  and  added  to  the  model's  crack  length 


a 


c 


+  CE  -  CFE 
Ac/Aa 


(57) 


where  c/a  is  determined  from  the  model  for  nearby  crack  lengths.  The 
compliance  Cp^  pertains  to  a  model  crack  length  of  aQ. 


The  effective  stress  intensity  factor,  «efp,  in  Table  4  was  calcu¬ 
lated  from  Equation  56  using  the  crack  length,  am,  and  the  load  given 
in  Table  3.  The  stress  intensity  factor,  KpE,  is  also  calculated  for 
convenience  from  Equation  56  but  using  the  crack  length  aQ.  For  another 
check  on  the  element  mesh,  the  elastic  stress  intensity  factor,  Kpp, 

was  also  calculated  based  on  J  integral  values  determined  from  VISCO 
for  several  paths.  For  linear  elastic  plane  stress  behavior  (Reference 
73) 

Kj  =  /EJ  (58) 

(Appendix  C  describes  the  VISCO  routine  for  calculating  the  J  integral). 


Figure  36a  shows  a  scaled  drawing  giving  four  different  J  integral 
paths  used  in  VISCO  for  the  center  cracked  plate  specimen.  Figure  36b 
shows  normalized  stress  intensity  factors  calculated  from  Equation  57 
and  the  J  values  along  paths  one  through  four  in  Figure  36a.  These  J 
values  were  from  a  linear  elastic  VISCO  analysis.  For  linear  elastic 
material,  the  J  integral  is  theoretically  path  independent  and  this  path 
independence  is  demonstrated  in  Figure  36.  It  should  be  noted  that  the 
good  agreement  between  VISCO  results  and  Equation  56  in  Figure  36b  also 
indicates  that  the  finite  element  mesh  (i.e.,  Figure  28)  accurately 
represents  the  center  cracked  plate. 
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Figure  36b.  Stress  Intensity  Factors  from  VISCO  J  Integrals  &  Equation  57 
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Table  5  summarizes  the  basic  details  of  the  VISCO  computer  runs 
which  employed  the  HEN  procedure.  The  computer  run  numbers  designated 
by  SI  through  S7  will  frequently  be  referred  to  in  the  following  discus¬ 
sions.  As  described  in  Section  IV,  each  of  these  runs  incorporated  the 
Bodner  material  model  and  the  crack  was  grown  in  the  model  at  a  sufficient 
rate  such  that  model  displacements  matched  test  data  with  time.  Each 
node  released  was  unloaded  in  five  seconds  except  runs  SI  and  S7  as 
indicated  in  Table  5.  Table  6  summarizes  the  basic  details  for  VISCO 
runs  similar  to  Table  5  but  with  no  crack  growth  allowed.  It  can  be 
seen  from  these  tables  that  the  computer  time  required  for  the  high  load 
runs  (e.g.,  S2  and  A3)  is  much  higher  than  for  the  lower  K  levels.  In 
the  case  of  run  A3,  computer  time  required  is  high  due  to  a  large  load 
causing  a  great  amount  of  plastic  flow  to  occur.  Recall  :  om  Section  III 
that  high  plastic  strain  rates  result  in  small  time  step  size  which  then 
requires  more  times  steps  to  simulate  a  given  amount  of  time  relative  to 
the  case  of  low  plastic  strain  rates.  Furthermore,  in  the  case  of  run 
S2,  extensive  crack  growth  occurred  requiring  19  nodes  to  be  released. 

Each  node  release  also  requires  relatively  small  time  steps  due  to  the 
redistribution  of  stress  and  the  associated  plastic  straining  around  the 
crack  tip. 

Figures  37  through  42  show  the  match  of  VISCO  displacements  with 
each  particular  test's  displacement  data.  Also,  the  amount  of  displace¬ 
ment  in  VISCO  for  no  crack  growth  under  load  is  given.  These  displace¬ 
ments  are  relative  to  the  displacements  existing  at  the  time  maximum  load 
is  achieved  which  means  all  displacements  prior  to  reaching  maximum  load 
were  not  included  in  this  test  data.  Maximum  load  was  normally  achieved 
in  five  seconds.  The  experimental  optical  technique  for  measuring 
displacements  is  highly  sensitive  and  can  resolve  displacements  in  the 
neighborhood  of  0.1  micron.  However,  this  technique  loses  sensitivity 
when  displacement  measurements  are  much  larger  than  a  micron.  Since  the 
main  interest  was  to  record  displacements  after  reaching  maximum  load 
(i.e.,  the  creep  crack  growth  data),  and  due  to  the  desire  to  maintain 
high  measurement  resolution,  the  load  application  displacements  were 
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Figure  38a.  Test  6  HEN  VISCO  Results,  Displacements  Matched  to  Test  Data 
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Figure  38b.  Test  6  HEN  VISCO  Results,  Crack  Growth 
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Figure  40a.  Test  9  HEN  VISCO  Results,  Displacements  Matched  to  Test  Data 


Figure  40b.  Test  9  HEN  VISCO  Results,  Crack  Growth 


AFWAL-TR-80-4140 


TIME  (min) 


Figure  41a.  Test  12a  HEN  VISCO  Results,  Displacements  Matched  to  Test 
Data 
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Figure  41b.  Test  12b  HEN  VISCO  Results,  Crack  Growth 
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Figure  42a.  Test  12b  HEN  VISCO  Results,  Displacements  Matched  to  Test 
Data 
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Figure  42b.  Test  12b  HEN  VISCO  Results,  Crack  Growth 
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left  out  of  the  creep  crack  growth  displacement  data.  Test  6  was  an 
exception  where  displacement  measurements  began  after  five  minutes  of 
test  time  had  expired. 

The  main  details  of  the  HEN  procedure  can  be  graphically  seen  In 
Figures  37a  and  b.  Observe  in  Figure  37a  the  two  specific  curves  of 
crack  tip  displacement.  One  from  VISCO  with  no  crack  growth  and  the 
other  being  test  data.  These  curves  start  to  deviate  from  each  other 
at  a  time  of  approximately  four  minutes  and  for  a  time  of  20  minutes 
the  difference  is  fairly  great  as  is  shown  by  a  bracket  in  Figure  37a. 
This  bracketed  difference  is  attributed  to  physical  crack  growth  and 
requires  release  of  crack  tip  nodes  in  this  simulation. 

Figures  37b  through  42b  depict  the  resulting  crack  growth  from  the 
HEN  VISCO  runs  for  each  test.  Note  that  due  to  this  discrete  finite 
element  analysis  technique,  the  crack  growth  is  a  step  function  (i.e.  , 
release  of  individual  crack  tip  nodes)  whereas  realistically  the  crack 
might  in  general  be  growing  in  a  smoother  manner  with  time. 

a.  Comparison  of  Results  from  Using  Different  Node  Unloading 
Methods 

In  cases  where  the  total  creep  crack  growth  is  only  a  few  node 
distances  (7.81  x  10”4  in. =20  microns)  the  displacements  developed  by 
VISCO  deviate  significantly  from  the  test  displacement  versus  time 
curve  (e.g.,  see  Figure  41a)  which  Implies  that  the  model  is  too  com¬ 
pliant  or  the  unloading  of  the  crack  tip  nodes  is  too  rapid.  These 
deviations  become  less  significant  for  larger  amounts  of  crack  growth 
as  seen  in  Figure  37a.  One  approach  to  minimize  these  deviations  for 
small  crack  growth  cases  is  to  make  a  second  VISCO  run  for  the  same  test. 
In  this  second  run,  node  release  times  from  the  first  VISCO  run  are 
Input.  Therefore,  the  continuous  node  unload  method  can  then  be  used 
as  described  in  Section  IV. 

Figure  37a  shows  VISCO  displacements  from  employing  the  five  second 
node  unload  method  (VISCO  Run  S3)  and  the  continuous  node  unload  method 
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(VISCO  Run  S7).  Note  that  the  displacements  from  Run  S7  are  only  slightly 
below  and  run  parallel  to  the  test  data  with  no  excursions.  Also,  Run 
SI  in  Figure  39a  matches  test  data  with  no  excursions  and  this  is  a  case 
similar  to  Figure  41a  where  overall  crack  growth  was  small  and  deviations 
were  large. 

Later  on  crack  growth  criteria  will  be  postulated  based  on  stress 
and/or  strain  at  the  crack  tip  and  therefore  reliable  data  from  the  HEN 
procedure  is  required  for  these  two  quantities.  Figures  43  and  44  show 
the  differences  in  stress  and  plastic  strain  at  the  crack  tip  which 
develop  between  VISCO  runs  using  the  two  different  node  unload  methods. 

For  large  amounts  of  crack  growth  in  this  case,  the  relative  difference 
in  effective  stress  at  the  crack  tip  is  seen  to  reach  a  steady  state 
value  of  approximately  7%.  The  relative  differences  of  the  y-component 
of  plastic  strain  at  the  crack  tip,  as  seen  in  Figure  44  is  much  less. 
These  differences  in  stress  and  strain  are  relatively  small  and  it  is 
concluded  that  essentially  the  same  results  ( i . e . ,  stress  and  strain 
at  the  crack  tip)  can  be  achieved  for  either  node  unload  method  for 
cases  with  large  amounts  of  crack  growth.  For  crack  growth  under  small 
loads  or  for  small  crack  growth  rates  more  time  is  allowed  for  stress 
relaxation  ahead  of  the  crack  tip.  Thus  on  the  average  crack  tip  stress 
and  strain  also  differ  little  between  the  two  unload  methods  for  slow 
crack  growth.  With  this  as  background,  it  was  decided  that  all  subse¬ 
quent  solutions  would  be  carried  out  using  the  five  second  node  unload¬ 
ing  scheme. 

b.  Dependence  on  Deformation  History 

Tests  8b  through  8d  were  not  included  in  the  HEN  VISCO  runs  since 
no  fatigue  precracking  was  done  prior  to  these  tests.  Without  fatigue 
precracking  these  tests  were  considered  to  be  a-typical  due  to  their 
different  prior  deformation  history.  To  check  out  dependence  on  prior 
deformation  history  the  following  VISCO  analysis  was  performed. 

Three  VISCO  runs  with  no  crack  growth  were  made  to  simulate  Test  8b 
and  its  dependence  on  prior  deformation.  Each  run  had  one  of  three 
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Figure  43.  Crack  Tip  Effective  Stress  at  the  Time  the  Crack  Tip  Node  Is 
Begun  to  be  Released 


AFWAL-TR-80-4140 


Crack  Tip  Plastic  Strain  at  the  Time  the  Crack  Tip  Node  Is 
Begun  to  be  Released 
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load-time  profiles  in  Figure  45  input  to  it.  The  VISCO  run  A7  had  no 
prior  load  hsitory,  run  A8  had  Test  8a  load  history  with  a  total  unload 
and  re-load  cycle,  which  is  the  most  realistic  load-time  model  of  Tests 
8a  and  b,  and  finally  run  A9  had  Test  8a  load  history  with  no  unload, 
but  a  load  increase  to  the  Test  8b  load  level.  The  top  of  Figure  45 
shows  the  displacement  versus  time  profile  for  each  of  these  VISCO  runs 
compared  to  Test  8b  data.  Again  these  displacements  are  relative  to  the 
displacements  existing  when  the  maximum  Test  8b  load  is  achieved. 


In  the  case  of  run  A8  where  complete  unloading  occurs  prior  to  Test 
8b  load,  the  VISCO  results  showed  that  plastic  flow  reoriented  itself 
during  the  unloading  due  to  the  reversing  of  the  principal  stress  at  the 
crack  tip  to  compression.  This  stress  develops  during  unloading  since 
the  material  has  prior  tensile  plastic  strains  from  the  Test  8a  load  and 
cannot  return  to  its  original  strain  free  state.  This  compressive 
behavior,  if  associated  with  crack  growth,  leads  to  the  crack  closure 
phenomenon  described  by  Newman  (Reference  83). 

Comparing  displacements  of  runs  A8  and  A9,  where  no  unloading  was 
done,  indicates  quite  similar  behavior.  When  no  prior  deformation  history 
exists  as  for  run  A7  the  displacements  differ  significantly,  as  shown  in 
Figure  45,  from  those  with  prior  load-deformation  history.  Note  that 
run  A7  fits  the  test  data  curve  best.  However,  test  compliance  data 
indicates  that  some  crack  growth  occurred  in  Test  8b.  If  displacements 
associated  with  this  crack  growth  were  included  in  these  VISCO  runs, 
experience  with  Test  8a  would  indicate  that  runs  A8  and  A9  would  be 
brought  up  to  the  test  data  and  run  A7  would  be  much  in  excess  of  test 
results.  It  should  also  be  noted  that  the  form  of  the  Bodner  model 
employed  in  VISCO  does  not  include  the  Bauschinger  effect  prevalent  in 
metals  when  compressive  yielding  follows  tensile  yielding.  Therefore, 
for  the  current  investigation  Tests  8b  through  8d  will  not  be  further 
studied  in  order  to  minimize  dependence  on  pretest  deformation  history 
in  this  investigation. 


100 


AFWAL-TR-80- 41 40 


oVISCO  Run  A 7 


vVISCO  Run  A8 


AVISCO  Run  A9 


Figure  45.  Load  History  Dependence  of  Displacements 
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c.  Crack  Growth  Results 

Table  7  displays  the  total  creep  crack  growth  increments  measured 
and  calculated  for  the  indicated  test.  The  column  entitled  Aam  presents 
the  measured  crack  growth  based  on  the  average  of  four  measurements 
across  the  thickness  as  shown  in  Figure  45.  The  column  entitled  Aa-| 
is  the  total  crack  growth  calculated  in  the  HEN  VISCO  runs  and  shown 
previously  in  Figures  37b  through  42b.  Two  other  convenient  methods  for 
calculating  crack  growth  without  the  need  of  HEN  VISCO  runs  were  also 
used  and  their  results  are  presented  as  Aa£  and  Aa3>  The  method  for 
calculating  Aa^  is  an  original  technique  developed  in  the  present 
investigation  and  will  be  discussed  subsequently. 

The  column  entitled  Aa-j  is  the  crack  growth  calculated  by  Clarke's 
elastic  compliance  method  (Reference  75).  This  method  predicts  negative 
crack  growth  in  some  instances  where  compliance  was  observed  to  decrease 
as  discussed  in  Section  IV.  Figure  46  compares  crack  growth  calculated 
by  Clarke's  method  to  the  results  from  a  HEN  VISCO  run  for  Test  9.  Note 
that  experimental  elastic  compliance  changes  from  Test  9  indicate  some 
unrealistic  negative  crack  growth  while  incorporating  Clarke's  method 
whereas  the  VISCO  results  show  a  realistic  monotonical ly  increasing 
amount  of  crack  growth  with  time.  However,  for  times  greater  than 
twenty  minutes,  both  curves  are  approximately  parallel  which  lends 
support  to  Clarke's  method  for  large  amounts  of  creep  crack  growth 
also  demonstrated  by  Donat  (Reference  76). 

The  column  in  Table  7  entitled  Aa^  is  the  crack  growth  calculated 
by  a  variation  of  Clarke's  elastic  compliance  method  as  described  below. 
In  this  variation  no  unload/reload  cycle  data  are  necessary  during  a 
creep  crack  growth  test  as  discussed  in  Section  IV.  The  VISCO  simula¬ 
tion  of  the  test  is  also  simplified  since  a  VISCO  run  is  made  using  test 
conditions,  but  no  crack  growth  is  allowed  in  the  VISCO  model.  Thus  any 
Increase  in  displacements  after  reacning  maximum  load  in  the  no-crack- 
growth  VISCO  run  can  only  occur  due  to  time  dependent  plastic  deformation 
allowed  by  the  Bodner  material  model.  Examples  of  these  no-crack-growth 
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iOnly  modeled  part  of  test  due  to  insufficient  data. 
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Figure  46.  Comparison  of  Crack  Growth  from  Experimental  Compliance  with 
HEN  VISCO  Results 
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VISCO  displacements  were  given  in  Figures  37a  through  42a.  Through 
simultaneous  use  of  test  displacement  data  and  the  no-crack- growth 
VISCO  results,  the  crack  growth  Aaj  Is  calculated  as  follows 


where  <$yest  is  the  experimental  displacement. 


The  displacement 


is  the  value  calculated  in  VISCO  for  the  same  test  conditions  but  no 


crack  growth  is  allowed.  The  denominator  term,  Ac/Aa  is  the  rate  of 
change  of  elastic  compliance  with  crack  length  generated  from  elastic 
VISCO  runs  of  the  test  specimens.  By  releasing  crack  tip  nodes  and 
dividing  the  resulting  VISCO  displacements  by  the  load  to  get  the 
respective  compliances,  the  curve  of  compliance  versus  crack  length 
was  determined  for  each  test's  indent  or  displacement  measurement 
location  (Figure  32).  The  slope  of  this  compliance  versus  crack 
length  curve  then  provided  Ac/Aa  ratio.  The  applied  load  P  is  divided 
into  the  difference  of  <5jest  and  6y  to  generate  a  change  in  compliance. 


This  change  in  compliance  is  assumed  to  be  primarily  dependent  on  crack 
growth  which  means  that  the  plastic  zone  developing  around  this  fixed 
crack  tip  is  approximately  the  same  size  for  an  extending  crack  and 
simply  translates  along  with  the  crack  tip. 


Figure  47  graphically  defines  6yest 


and  <5 

v 


at  an  arbitrary  measure¬ 


ment  time  t  .  Note  that  this  Aa,  method  is  effectively  the  same  as  the 
me 


HEN  VISCO  procedure  for  small  amounts  of  crack  growth  where,  as  mentioned 
previously,  the  plastic  zone  size  is  assumed  constant.  Therefore  a 
curve  of  crack  growth  versus  time  could  also  be  generated  by  this  Aa^ 
method  by  applying  Equation  59  continuously  along  the  test  data  curve. 


Consider  the  measured  crack  growth  column  Aa^  and  the  HEN  VISCO 
results  Aa-| .  The  differences  between  these  two  columns  is  only  +  10% 
for  those  tests  totally  modeled.  Only  part  of  Test  6  was  modeled  since 
the  number  of  uniform  crack  growth  elements  ahead  of  the  initial  crack 
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Figure  47.  Schematic  of  Displacements  Used  to  Calculate  Crack  Growth  by 
Equat’"'  “  59 
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tip  were  already  exhausted  after  half  of  the  test  duration.  This  also 
occurred  in  modeling  Test  9.  Hence,  in  both  Tests  6  and  9  the  value  of 
Aa-j  was  determined  by  extrapolating  the  final  or  relatively  steady  state 
slope  of  the  crack  growth  versus  time  curves  in  Figures  38b  and  40b. 
Another  special  characteristic  of  Test  6  is  that  displacement  data  versus 
time  was  not  measured  during  the  first  five  minutes  of  the  test.  Thus, 
based  on  other  test  behavior,  a  significant  amount  of  the  early  relatively 
rapid  displacement  rates  were  ignored.  Some  of  this  early  displacement 
would  apparently  have  required  more  crack  growth  in  VISCO  than  the  exist¬ 
ing  data  did  while  using  the  HEN  procedure. 

The  crack  growth  values,  Aa2>  in  Table  7  also  correlate  quite  well 
with  the  measured  values  Aam.  The  good  correlation  of  Aa2  with  test  data 
provides  support  to  this  convenient  method  developed  herein  for  calculat¬ 
ing  creep  crack  growth.  The  significant  convenience  feature  in  calculat¬ 
ing  Aa2  is  that  once  a  no-crack-growth  VISCO  run  is  made  for  a  given  set 
of  test  conditions.  Equation  59  can  simply  be  applied  for  analyzing  all 
further  experiments  with  approximately  the  same  test  conditions  (e.g., 
load,  crack  length  and  geometry). 

A  somewhat  similar  HEN  analysis  was  done  by  HSU  et.  al  (Reference 
32)  to  predict  crack  growth  in  zirconium  at  elevated  temperatures.  HSU's 
crack  growth  predictions  were  2.5  times  greater  than  actual  test  data. 

In  the  present  HEN  analysis,  crack  growth  predictions  were  within 
approximately  10%  of  test  data  for  the  tests  that  were  totally  simulated. 
The  extremely  good  correlation  using  the  present  method  is  attributed  to 
both  higher  resolution  experimental  displacement  data  than  HSU's  and  a 
more  realistic  material  model  that  includes  creep  behavior.  Additional 
displacements  due  to  creep  in  HSU’s  analysis  would  have  reduced  the 
amount  of  crack  growth  predicted.  Since  HSU’s  predictions  were  high, 
this  reduction  would  be  in  the  direction  of  better  correlation  with 
test  data. 

Another  comparison  of  current  results  can  qualitatively  be  made  with 
Newman’s  finite  element  analysis  of  fatigue  crack  propagation  (Reference 
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83).  Figure  48  shows  crack  opening  displacement  profiles  during  creep 
crack  growth  from  the  VISCO  model  as  represented  by  the  solid  lines. 

Note  that  with  the  crack  tip  at  point  zero,  the  elastic  crack  opening 
displacement  profile  (f.e.,  the  dashed  line  In  Figure  48)  is  lower  than 
the  solid  line  which  incorporates  the  Bodner  material  model.  However, 
as  the  crack  grows  a  wake  of  residual  plastic  deformation  Is  left  behind 
the  crack  tip.  After  12  increments  of  crack  growth.  Figure  48  shows 
how  this  plastic  wake  diminishes  the  elastic-plastic  crack  profile  below 
a  purely  elastic  profile  for  the  same  crack  length  "a".  The  residual 
plastic  deformation  indicated  by  the  cross-hatched  area  in  Figure  48 
was  also  displayed  in  a  similar  figure  by  Newman. 

2.  CRACK  GROWTH  CRITERION 

Based  on  the  good  correlation  between  actual  and  predicted  crack 
growth,  the  VISCO  results  from  the  HEN  applications  were  examined  for 
potential  crack  growth  parameters.  The  local  crack  tip  parameters  such 
as  strain  and  crack  opening  displacement  (C.O.D.)  were  examined  initially 
since  these  parameters  have  shown  promise  elsewhere  (References  84,85) 
for  correlating  finite  element  results  with  crack  growth  test  data. 

The  main  goals  here  were  to  check  out  the  validity  of  existing  crack 
growth  criteria  (i.e.,  critical  C.O.D.  and  strain),  possibly  modify  one 
of  these  existing  criteria  to  better  fit  test  data,  and  postulate  a  new 
criterion  that  might  better  account  for  creep  damage  accumulation  and 
crack  growth  displayed  by  the  HEN  VISCO  analysis. 

a.  Critical  Strain  Criterion 

Examination  of  the  strains  in  the  elements  adjacent  to  the  crack 
tip  prior  to  each  node  release  for  crack  growth  in  the  HEN  VISCO  runs 
revealed  that  no  single  value  for  the  critical  strain  would  satisfy  all 
test  conditions.  However,  a  few  VISCO  runs  were  made  using  a  fixed 
critical  strain  criterion  to  further  evaluate  its  applicability. 
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Figure  48.  Crack  Opening  Displacement  Profiles  During  Creep  Crack  Growth 
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The  critical  strain  criterion  was  Implemented  in  VISCO  by  comparing 
the  average  of  the  plastic  strain  components  normal  to  the  crack  and 
within  elements  adjacent  to  the  crack  tip  with  a  critical  strain  value, 
ecrit’  as  *1me  Pr09resse<*-  When  the  average  crack  tip  plastic  strain 

exceeded  the  critical  value,  ecr1t»  the  crack  tip  node  was  released  and 

unloaded  in  five  seconds.  Table  8  gives  a  summary  of  the  basic  details 
for  VISCO  runs  employing  the  critical  strain  criterion.  A  ecr^t  value 

of  0.030  was  found  to  work  well  from  the  HEN  VISCO  runs  of  Test  8a  where 
KpE  was  16.3  ksi/TrT.  However  for  Test  9  where  KpE  was  36.8  ksi/Tn  the 

e  value  needed  to  be  0.090  to  work  well  as  shown  in  Figure  49a.  The 

corresponding  ecr^t  dependent  crack  growth  versus  time  is  given  in  Figure 

49b.  The  ecrit  value  of  0.075  allowed  too  much  displacement  compared 

to  test  data  as  seen  in  Figure  49.  The  displacement  results  for  a  e 

value  of  0.090  appear  to  fit  the  data  quite  well  over  the  first  five 
minutes  of  the  test.  Likewise  the  resulting  crack  growth  in  Figure  49b 
agrees  quite  well  with  the  HEN  data  in  Figure  39b.  However  based  on 
examination  of  the  prior  HEN  VISCO  runs,  if  the  simulation  time  were 
continued,  the  ecrit  value  of  0.090  would  have  been  too  large.  Therefore, 

Insufficient  crack  growth  would  have  resulted  and  the  displacements 
would  have  fallen  away  from  the  test  data  as  shown  by  the  extrapolated 
dashed  curve. 

Figure  50  displays  the  HEN  VISCO  crack  tip  strain  values  taken  at 
the  time  a  crack  tip  node  was  to  be  released  in  the  HEN  VISCO  runs. 

This  plot  was  motivated  in  an  attempt  to  find  a  general  trend  of  the 
critical  strain  values  or  to  determine  a  mathematical  relationship  with 
time  to  envelope  the  HEN  results.  Note  that  in  general  the  strains  are 
high  for  short  times  and  then  diminish  with  time  to  a  common  value  of 
approximately  0.03.  In  order  to  develop  a  critical  strain  functional 
expression  one  can  see  a  need  for  a  decaying  parameter  which  would  fit 
the  upper  bound  strain  values  in  Figure  50  related  to  rapid  crack 
growth.  In  addition  the  expression  must  have  the  capability  of  allowing 
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TIME  (min) 


Figure  49a.  VISCO  Critical  Strain  Crack  Growth  Criterion,  Resulting 
Displacements  Due  to  Crack  Growth 


0  5  10  15 


TIME  (min) 

Figure  49b.  VISCO  Critical  Strain  Crack  Growth  Criterion,  Resulting 
Crack  Growth 
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Figure  50.  Time  Dependent  Critical  Strain  Crack  Growth  Criterion 


113 


AFWAL-TR-80-4140 


the  region  around  the  crack  tip  node  to  deteriorate  with  exposure  to  the 
crack  tip  environment  (i.e.,  not  only  the  environment  external  to  the 
specimen  but  also  the  singular  strain  state  internal  to  the  specimen). 
Crack  tip  deterioration  is  considered  to  be  displayed  by  the  lower  bound 
of  strain  values  in  Figure  50.  Therefore,  the  general  properties  of  the 
critical  strain  function  have  been  stated.  The  curve  must  be  decaying 
with  test  time,  thus  a  negative  exponential  function  is  in  order.  And 
since  critical  strain  values  appear  to  diminish  with  crack  tip  exposure 
time,  but  not  as  rapidly  as  an  exponential  function  would  dictate,  the 
cosine  function  was  examined. 


An  empirical  expression  for  the  critical  strain  that  fits  the  HEN 
VISCO  results  fairly  well  is  given  by 


|  eQ[A  exp(-bt)  cos/|y  V  t 
"crit  I  0 


if  T  <  T 
—  o 


if  T  >  T 


(60) 


-1 


where  e  =  0.03,  A=3.0,  b=l .34x10  sec  ,  T  =600. sec. 
o  o 


The  value  of  eQ  represents  the  critical  strain  for  large  time.  The  co¬ 
efficient  A  is  determined  at  t  equal  to  zero.  Once  A  and  e  are  chosen 

o 

b  is  determined  by  best  fitting  the  upper  bound  where  crack  tip  exposure 
time  is  small  and  set  to  zero  (i.e.,  T  =  o).  The  parameter  Tq  determines 

how  rapidly  the  critical  strain  diminishes  to  eQ  with  crack  tip  exposure 

time  T. 


Motivation  for  the  development  of  Equation  60  comes  from  environ¬ 
mental  effects  such  as  oxidation  at  these  high  temperatures  from  exposure 
to  laboratory  air.  This  oxidation  is  then  associated  with  changing 
material  properties  at  the  crack  tip  such  as  the  critical  strain  value. 
Crack  growth  in  alloys  similar  to  IN-100  has  been  found  to  be  quite 
sensitive  to  the  environment  and  the  resulting  oxidation  that  can  occur 
when  an  elevated  temperature  crack  growth  test  is  done  in  air  (Reference 
86).  Therefore,  Equation  60  is  an  attempt  to  represent  the  rate  at  which 
the  critical  strain  for  crack  growth  is  diminished  with  time  due  to 
environmental  effects. 
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The  critical  strain  criterion  just  formulated  would  be  employed  in 
VISCO  in  the  following  steps. 

1.  Register  both  total  test  simulation  time,  t,  and  crack  tip 
exposure  time,  T. 

2.  Evaluate  the  critical  strain  level  during  each  time  step  from 
Equation  60. 

3.  Compare  with  crack  tip  element's  plastic  strain  accumulating  in 
the  VISCO  analysis. 

4.  Release  crack  tip  node  when  VISCO  plastic  strain  at  the  crack 
tip  exceeds  e  from  Equation  60. 

5.  Crack  tip  exposure  time  T  is  set  to  zero  and  the  above  steps  1-4 
are  repeated  for  the  next  crack  tip  node. 

b.  Critical  Crack  Opening  Displacement  Criterion 

The  crack  opening  displacement  is  defined  here  as  the  C.O.D.  at  the 
first  node  behind  the  crack  tip  in  the  finite  element  model.  Examination 
of  C.O.D. 1 s  taken  from  HEN  VISCO  runs  in  Table  5  just  prior  to  releasing 
a  crack  tip  node  revealed  no  single  C.O.D.  value  that  could  be  used  for 
all  test  conditions  as  a  critical  C.O.D.  A  value  of  0.280  x  10“4  inches 
was  found  from  HEN  VISCO  results  for  Test  8a  with  KpE  =  16.3  ksi/fn  to 

work  best,  yet  when  a  few  VISCO  runs,  as  summarized  in  Table  9,  were  done 
with  a  critical  C.O.D.  criterion  and  an  increase  to  KfE  =  36.8  ksi/irT 

for  Test  9,  the  best  critical  C.O.D.  became  approximately  0.500  x  10"^ 
inches  as  shown  in  Figure  51.  Further  observance  of  Table  9  and  Figure 
51  shows  that  this  criterion  was  very  sensitive  to  small  changes  in 
critical  C.O.D.  Therefore  further  evaluation  of  this  criterion  was 
deemed  unnecessary  as  compared  with  the  less  sensitive  critical  strain 
criterion. 

c.  Critical  Damage  Accumulation  Criterion 

In  several  theoretical  works  and  as  stated  recently  by  Goodall  and 
Chubb  (Reference  62),  creep  rupture  of  uncracked  components  under  a 
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Figure  51.  VISCO  Critical  C.O.D.  Criterion,  Resulting  Displacements 
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varying  stress  history  is  governed  by  the  life  fraction  rule.  Conse¬ 
quently,  this  reference  indicated  that  for  a  uniaxial  stress  history, 
o(t),  rupture  occurs  at  a  time  t  given  by 


/k_  -  , 

0  tr(o) 


(61) 


where  tr(o)  is  the  rupture  time  corresponding  to  a  constant  stress  level, 

a.  When  experimental  stress  values  are  plotted  against  their  rupture 
times  on  logarithmic  scales  the  relationship  is  often  linear  in  the 
region  of  practical  interest.  Thus  if  M  and  C  are  material  constants 
it  is  assumed  that 

oM  tr(o)  =  C  (62) 

Substituting  Equation  62  into  61  yields 

rlr 

J  oM  dt  =  C  (63) 

o 


The  description  of  creep  rupture  given  by  Equations  61  through  63  was 
discussed  by  Goodall  and  Chubb.  The  reference  mentions  that  this  rupture 
model  is  only  one  of  several  possible  formulations.  However,  neither 
experimental  nor  theoretical  work  has  provided  an  alternative  to  Equation 
61  that  gives  a  significantly  better  description  of  material  response. 


It  is  recognized  that  Equation  62  applies  most  directly  to  creep 
rupture  of  uncracked  uniaxial  components.  However,  the  present  author 
considers  that  similar  behavior  might  be  possible  in  creep  crack 
propagation.  A  schematic  of  the  postulated  behavior  involved  with  creep 
crack  propagation  is  given  in  Figure  52.  This  figure  shows  a  creep 
damage  front  preceding  the  crack.  Within  this  front,  the  material  is 
accumulating  creep  damage  in  the  form  of  microcracks.  This  type  of 
creep  damage  is  also  associated  with  creep  rupture  of  uncracked 
components. 
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I 

!  Figure  52.  Schematic  Representation  of  Creep  Crack  Propagation 
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In  order  to  apply  Equation  63  as  a  crack  growth  criterion,  a  process 
zone  6  is  required  and  defined  in  Figure  52.  In  addition,  the  rupture 
time  t  is  redefined  from  Equation  61  to  the  elapsed  time  the  crack 
requires  to  grow  from  one  node,  in  the  HEN  VISCO  results,  to  the  next 
node.  In  other  words,  it  is  the  time  period  during  which  the  process 
zone  6  Is  exposed  to  the  crack  tip  stress  field  prior  to  rupture.  In 
the  VISCO  finite  element  analysis,  this  process  zone  was  taken  as  one 
element  preceding  the  crack  tip.  The  average  component  of  stress 
normal  to  the  crackline  from  three  elements  adjacent  to  the  crack  tip 
was  used  as  the  stress,  <$,  in  Equation  63  as  shown  in  Figure  52. 


Since  the  greatest  stress  exists  at  the  crack  tip  and  environmental 
degradation  is  considered  to  be  most  prevalent  there,  most  damage  accumu¬ 
lation  was  assumed  to  occur  in  the  process  zone  after  the  arrival  of  the 
crack  tip  to  the  process  zone's  border.  Therefore  time  in  Equation  63 
was  measured  from  crack  tip  arrival  time  to  the  current  crack  tip  node. 


or 


/•vs 

j  oM 

*A 


dt  =  C 


(64) 


The  constants  M  and  C  were  determined  based  on  results  from  the  HEN  VISCO 
runs.  To  accomplish  this.  Equation  64  was  approximated  as 


(65) 


The  rupture  time  or  crack  growth  times  t  were  taken  from  HEN  VISCO 
results.  The  stress,  cgVg,  also  based  on  HEN  results,  was  an  average 

over  time  t  of  the  crack  tip  stress  defined  previously.  Since  the 
interest  here  is  to  develop  values  for  M  and  C  which  apply  to  the 
entire  set  of  tests,  it  became  obvious  that  there  were  more  combinations 
of  tr  and  <?aVg,  than  necessary  to  uniquely  define  M  and  C.  Consequently, 

to  include  the  data  from  each  HEN  VISCO  run,  a  least  square  fit  of  the 
data  on  a  log-log  plot  was  used  to  best  fit  the  time-stress  data.  The 
values  chosen  were 

C=8.63  x  1079  (psi)^  sec 
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Equation  64  along  with  the  previous  constants  were  then  incorporated 
into  VISCO  as  a  critical  damage  accumulation  criterion  for  crack  growth. 
Table  10  summarizes  the  basic  details  of  the  VISCO  runs  using  this  crack 
growth  criterion. 

Figures  53  through  58  give  the  results  of  applying  VISCO  with  the 
critical  damage  accumulation  criterion  to  the  indicated  test  number  in 
each  figure.  The  resulting  VISCO  displacements  are  compared  to  test 
data  in  Figures  53a  through  58a  whereas  the  resulting  crack  growth  from 
VISCO  is  given  in  Figures  53b  through  58b.  For  the  higher  loadings  such 
as  Test  9  in  Figure  56a  the  difference  between  test  data  and  VISCO  results 
is  greatest.  Part  of  this  difference  may  be  due  to  the  fact  that  at  this 
high  load  (Kp^=36.8  ksi/Tn)  the  damage  zone  as  indicated  in  Figure  52  is 

larger  than  the  process  zone  used  herein  (i.e.,  one  element  size  or  a 
characteristic  dimension  of  7.81  x  10"^  inches).  Hence  significant 
damage  accumulation  may  occur  in  the  material  before  arrival  of  the 
crack  tip  or  time  t^  for  a  given  element  in  the  crack  path.  This 
damage  occurring  in  an  element  prior  to  tft  for  that  respective  element 
was  neglected  in  the  present  damage  accumulation  criterion.  However, 
for  the  lower  loads  agreement  with  the  test  data  is  quite  good  consider¬ 
ing  that  creep  rates  for  the  same  test  conditions  can  easily  vary  by  a 
factor  of  two  or  three  which  is  why  log-log  plots  are  used  to  plot  creep 
data  (Reference  29). 

It  should  also  be  noted  that  the  M  and  C  values  used  were  determined 
from  the  approximate  Equation  65  expression.  Crack  growth  results  from 
this  criterion  might  be  improved  by  iterating  or  making  small  changes 
to  M  and  C  and  making  further  VISCO  runs  in  an  effort  to  better  fit 
test  data. 

3.  CRACK  GROWTH  RATE  CRITERIA 

In  this  section  crack  growth  rate  criteria  will  be  discussed  based 
on  the  steady  state  crack  growth  rates  developed  by  the  HEN  VISCO  runs. 
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Figure  53a.  Test  5  -  VISCO  Critical  Damage  Accumulation  Criterion 
Displacement  Results 
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Figure  53b.  Test  5  -VI SCO  Critical  Damage  Accumulation  Criterion  Crack 
Growth  Resul ts 
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re  54a.  Test  6  -  VIS CO  Critical  Damage  Accumulation  Criterion 
Displacement  Results 


TIME  (min) 

Figure  54b.  Test  6  -  VISCO  Critical  Damage  Accumulation  Criterion  Crack 
Growth  Results 
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Figure  55a.  Test  8a  -  VISCO  Critical  Damage  Accumulation  Criterion 
Displacement  Results 
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Figure  55b.  Test  8a  -  VISCO  Critical  Damage  Accumulation  Criterion 
Crack  Growth  Results 
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Figure  56a.  Test  9  -  VISCO  Critical  Damage  Accumulation  Criterion 
Displacement  Results 
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Figure  56b.  Test  9  -  VISCO  Critical  Damage  Accumulation  Criterion  Crack 
Growth  Results 
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Figure  57a. 
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Figure  57b.  Test  12a  -  VISCO  Ci 
Crack  Growth  Resul 
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Figure  58a.  Test  12b  -  VISCO  Critical  Damage  Accumulation  Criterion 
Displacement  Results 


TIME  (min) 


Figure  58b.  Test  12b  -  VISCO  Critical  Damage  Accumulation  Criterion 
Crack  Growth  Results 
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The  steady  state  crack  growth  rates  were  determined  by  best  fitting  a 
straight  line  to  the  crack  growth  curves  In  Figures  37b  through  42b. 

The  Initial  transient  portion  of  these  crack  growth  curves  was  Ignored 
for  these  steady  state  values.  Also  a  discussion  of  the  merits  of  the 
*C  Integral  as  a  crack  growth  rate  criterion  will  be  given. 

a.  Stress  Intensity  Factor  Criterion 

Crack  growth  rate,  a,  has  previously  been  found  to  correlate  with 

the  elastic  stress  Intensity  factor  as  given  in  Equation  2  and  recalled 

here  .  a 

a  =  A  (K)  (66) 

This  relationship  plots  as  a  straight  line  on  log-log  paper,  as  shown  in 
Figure  59.  The  experimental  data  referred  to  in  Figure  59  was  for  IN-100 
behavior  at  1350°F  which  is  the  same  alloy  and  temperature  used  in  the 
present  investigation.  Note  that  Donat's  experimental  data  (Reference 
76)  covered  a  range  in  K  values  from  30  to  approximately  80  ksi  inches. 

In  order  to  compare  with  the  lower  K  levels  in  the  current  investigation, 
the  line  representing  the  best  fit  to  Donat's  data  was  extrapolated  as 
shown  by  the  dashed  line  in  Figure  59.  Agreement  with  the  present  HEN 
VISCO  results,  in  which  is  taken  from  Table  5  and  a  from  Figures 

37b  through  42b  is  good  especially  considering  the  fact  that  the  test 
data  line  was  extrapolated. 

This  criterion  has  the  distinct  advantage  relative  to  criteria  pre¬ 
sented  earlier  that  once  the  constants  A  and  a  are  determined,  it  can  be 
used  independent  of  finite  element  analyses.  This  advantage  is  due  to 
the  fact  that  K  can  be  calculated  for  most  test  geometries  by  relatively 
simple  equations  like  Equation  56.  Thus  the  so-called  steady  state 
creep  crack  growth  rate  can  be  simply  calculated  from  Equation  66,  but 
it  should  be  kept  in  mind  that  incubation  time  for  crack  initiation  nor 
the  initial  rapid  crack  growth  observed  in  the  HEN  VISCO  results  is 
captured  with  this  criterion. 
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Figure  59.  Crack  Growth  Rate  vs.  Stress  Intensity  Factor 
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b.  Net  Section  Stress  Criterion 

Crack  growth  rates  have  also  been  shown  to  be  related  to  net  section 
stress  as  given  in  Equation  3  and  recalled  here 

a  3  B(on)6  (67) 

Figure  60  shows  how  Equation  66  also  plots  as  a  straight  line  (the  dashed 
line)  on  log-log  paper.  The  present  results  seem  to  correlate  to  Equa¬ 
tion  67  as  well  as  they  did  to  Equation  65  on  log-log  paper.  This  is 
not  surprising  since  for  the  center  cracked  plate  it  can  be  shown  that 
K  is  approximately  directly  proportional  to  an  for  crack  lengths  of 
a/W  from  .2  to  .7  which  spans  the  crack  lengths  in  the  current  study. 

The  net  section  stress  criterion,  like  the  K  criterion,  also  neglects 
crack  growth  characteristics. 

c.  C*  Integral  Criterion 

The  C*  integral  is  an  extension  of  the  J  integral  concept  for 
application  to  creep  crack  growth  (a  detailed  discussion  of  C*  is  given 
in  Appendix  C).  Theoretically  the  C*  integral  is  path  independent  for 
a  creeping  solid  where  stress  is  only  a  function  of  the  plastic  strain 
rate  and  elastic  strain  rates  do  not  enter  into  the  formulation.  In  the 
present  research,  C*  was  calculated  as  shown  in  Appendix  C.  In  this 
calculation  total  strain  rates  were  attributed  to  creeping  plastic 
strain  rates  which  implies  no  elastic  strain  rates  exist.  For  a 
realistic  material  that  includes  elastic  behavior,  if  the  elastic  strain 
rates  are  zero  then  the  stress  state  is  constant  with  time,  and  there¬ 
fore,  the  creeping  plastic  strains,  since  they  are  a  function  of  stress 
must  also  be  constant.  Accordingly,  with  constant  stress  and  strain 
rate  the  W*  integral  was  integrated  directly  in  Appendix  C. 

The  previous  description  for  calculating  C*  in  VISCO  was  implemented 
and  applied  to  a  realistic  elastic-viscoplastic  material  (i.e.,  IN-100) 
being  studied  herein.  It  was  observed  that  during  load  application  the 
C*  values  were  extremely  high  due  to  the  elastic  strain  rate  contribution. 
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Figure  60.  Crack  Growth  Rate  vs.  Net  Section  Stress 
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After  maximum  load  was  achieved  C*  values  reduced  down  to  much  lower 
steady  state  values  until  crack  growth  began  in  the  HEN  VISCO  runs. 

Again  as  crack  growth  began,  the  C*  values  increased  significantly 
due  to  the  elastic  strain  rate  contribution  as  stress  was  redistributing 
around  the  moving  crack  tip.  It  should  be  noted  that  any  attempt  to 
remove  elastic  contributions  to  strain  and  displacements  will  result 
tn  an  ill -posed  problem  since  the  displacement  rate  components  needed 
In  Equation  7  cannot  be  resolved  Into  an  elastic  and  plastic  portion. 

Thus  far  evaluation  of  the  C*  integral  given  in  Equation  7,  the 
elastic  part  of  the  strain  rate  is  neglected,  otherwise  the  Integral 
has  no  meaning.  But  when  the  elastic  strain  rates  are  ignored,  C* 
should  only  be  calculated  after  the  stresses  are  fairly  constant  and 
prior  to  crack  Initiation,  C*  will  then  be  a  constant  until  the  crack 
begins  to  grow.  Therefore,  it  is  impossible  to  relate  C*  to  any 
incubation  time  for  crack  growth.  Moreover,  during  crack  growth  the 
contribution  from  elastic  strain  rates  is  again  substantial  and  cannot 
be  Ignored.  Thus,  It  appears  that  the  C*  integral  is  ineffective  as  a 
fracture  criterion  in  a  finite  element  model  for  creep  crack  growth 
in  the  current  investigation. 

d.  Load  Point  Displacement  Rate  Criterion 

An  expression  relating  crack  growth  rates  to  load  point  displace¬ 
ment  rate  was  given  in  Equation  4.  Unfortunately  this  criterion  suffers 
from  problems  similar  to  the  C*  integral.  The  load  point  displacement 
rate  after  reaching  maximum  load  is  the  sum  of  the  displacement  rate 
due  to  crack  growth  as  well  as  the  displacement  rate  due  to  plastic 
deformation.  Hence  the  load  point  displacement  rate  may  have  several 
values  for  the  same  crack  growth  rate  depending  on  the  rate  of  plastic 
deformation  such  as  demonstrated  in  Equation  59.  One  can  get  an  appre¬ 
ciation  of  how  much  plastic  deformation  contributes  to  the  overall 
displacements  by  looking  at  the  no-crack-growth  displacements  versus 
time  in  Figures  37a  through  42a  which  are  totally  a  result  of  plastic 
deformation.  Therefore  the  load  displacement  rate  does  not  provide  a 
unique  solution  to  the  crack  growth  rate  unless  variations  in  the  plas¬ 
tic  deformation  rate  can  be  neglected. 
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SECTION  VI 

SUMMARY  AND  CONCLUSIONS 

A  two-dimensional  (plane  stress/plane  strain)  finite  element  program 
has  been  developed  which  accounts  for  both  nonlinear  viscoplastic  material 
behavior  and  changing  boundary  conditions  due  to  crack  growth.  Three 
viscoplastic  material  models:  (1)  Malvern  Flow  Law,  (2)  Norton's  Creep 
Law,  and  (3)  Bodner-Partom  Flow  Law  were  incorporated  into  the  program. 
These  time  dependent  material  models  were  numerically  integrated  through 
time  by  a  linear  Euler  extrapolation  technique.  A  variable  time  step 
algorithm  was  included  that  maximized  time  step  size  during  the  analysis 
while  maintaining  good  accuracy.  This  program  was  used  as  the  plane 
stress  theoretical  model  for  the  hybrid  experimental -numerical  procedure 
employed  to  analyze  sustained  load  creep  crack  growth  test  data.  The 
test  specimens  were  center  cracked  plates  made  of  IN-100  and  tested  at 
1 350°F. 

The  following  statements  and  conclusions  are  based  on  the  creep 
crack  growth  analysis  herein. 

1.  A  method  for  getting  crack  growth  behavior  solely  from  displacement 
measurements  in  conjunction  with  a  cracked  specimen  model  which 
utilizes  realistic  constitutive  relationships  has  been  developed. 

The  constitutive  law  was  especially  tailored  to  the  nickel-base 
alloy  studied  which  displays  time  dependent  nonlinear  inelastic 
behavior  at  elevated  temperatures.  It  has  been  demonstrated  that 
the  technique  can  be  applied  where  crack  extension  is  very  small 
and  could  not  otherwise  be  resolved  by  conventional  experimental 
crack  measuring  techniques  (e.g.,  compliance  techniques  or  using  a 
travelling  microscope).  This  method  provides  realistic  mono- 
tonically  increasing  crack  growth  values  with  resolution  better 
than  0.001  inch.  Crack  growth  predictions  agreed  to  within  10% 
of  post-test  measurements. 
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2.  Numerical  procedures  were  developed  to  efficiently  Integrate  the 
nonlinear  time  dependent  material  models  and  simulate  crack  growth 
by  two  crack  tip  node  unloading  methods.  The  numerical  time  inte¬ 
gration  procedure  utilized  an  Euler  linear  extrapolation  technique 
with  a  variable  time  step  algorithm  that  maximized  time  step  size 
during  the  simulation  while  maintaining  good  accuracy.  One  crack 
tip  node  unload  method  diminished  the  restraining  force  on  the 
finite  element  crack  tip  node  in  a  specified  five-second  time 
period  independent  of  crack  growth  rate  whereas  the  second  method 
continuously  unloaded  crack  tip  nodes  in  proportion  to  a  predeter¬ 
mined  crack  growth  rate.  Unloading  the  nodes  continuously  provided 
a  closer  fit  to  displacement  versus  time  test  data  however  the 
average  displacement  versus  time  was  approximately  the  same  for 
both  node  unloading  methods. 

3.  A  procedure  was  developed  for  determining  crack  extension  using 
calculations  of  viscoplastic  deformation  with  no  crack  growth. 

In  this  procedure,  the  difference  between  total  test  deformation 
and  viscoplastic  deformation  is  attributed  to  crack  extension. 
Extremely  good  crack  growth  predictions  were  made. 

4.  The  elastic  compliance  method  for  resolving  creep  crack  extension 
has  been  shown  to  imply  negative  unrealistic  crack  growth  and  is 
unreliable  especially  during  the  first  part  of  a  creep  crack  growth 
test. 

5.  Several  parameters  were  studied  for  their  potential  as  creep  crack 
growth  controlling  parameters. 

a.  No  single  fixed  value  of  strain  for  a  critical  strain  crack 
growth  criterion  was  found  to  match  all  test  conditions  in  this 
investigation.  Environmental  effects  apparently  tend  to  lower 
the  critical  strain  magnitude  with  time,  under  load.  An 
empirical  relationship  was  developed,  based  on  the  HEN  results, 
which  gives  the  critical  strain  a  diminishing  value  with  time. 
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This  investigation  did  not  include  any  plane  strain  analyses. 
Due  to  higher  constraint  at  the  crack  tip  for  plane  strain, 
less  stress  redistribution  would  occur.  Therefore,  it  seems 
possible  that  the  critical  strain  for  all  test  conditions  may 
vary  less  in  a  plane  strain  simulation. 

b.  No  single  fixed  value  for  C.O.D.  was  found  to  match  all  test 
conditions  using  a  critical  C.O.D.  crack  growth  criterion. 

The  C.O.D.  behaved  similar  to  the  crack  tip  strain  with  time, 
however,  its  percent  variation  was  less. 

c.  A  critical  damage  accumulation  criterion  for  crack  growth  was 
developed  based  on  a  modification  of  the  life  fraction  rule 
for  creep  rupture  to  account  for  environmental  effects  at  the 
crack  tip.  Application  of  this  criterion  provided  good 
agreement  with  the  low  to  medium  load  test  conditions.  For 
the  highest  load  test  cases,  this  criterion  predicted  crack 
growth  rates  somewhat  lower  than  the  HEN  results.  It  appears 
that  accumulation  of  damage  over  all  time  and  not  just  crack 
tip  exposure  time  might  improve  the  results. 

6.  Data  obtained  in  this  investigation  through  numerical  calculations 
provided  crack  growth  versus  time  and  thus  crack  growth  rate  time, 
a.  The  a  data  compared  well  with  published  data  for  the  same 
material  and  temperature  when  plotted  against  stress  intensity 
factor.  The  present  data  was  obtained  for  a  and  K  values  lower 
than  the  referenced  data.  Net  section  stress  also  provided  good 
correlation  with  the  predicted  crack  growth  rates. 

7.  The  C*  integral  and  load  line  displacement  rate  were  investigated 
as  possible  parameters  controlling  crack  growth  rate,  a.  The  C* 
integral  is  an  unreliable  parameter  for  predicting  creep  crack 
growth  due  to  its  formulation  which  is  based  on  a  creeping  solid 
behavior  that  neglects  elastic  strain  rates.  The  load  line  dis¬ 
placement  rate  which  can  be  shown  to  be  proportional  to  C*  also 
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does  not  provide  a  unique  solution  for  the  crack  growth  rate  unless 
variations  in  plastic  deformation  rate  can  be  ignored.  In  general 
these  parameters  appear  to  have  no  applicability  to  crack  growth 
rate  prediction  using  numerical  modeling  of  materials.  However 
these  parameters  seem  to  correlate  a  data  fairly  well  once  the 
solution  is  known  as  seen  in  the  literature. 

8.  Constant  strain  triangular  finite  elements  of  the  size  of  grains 
at  the  crack  tip  work  well  for  resolving  small  increments  of  creep 
crack  growth  through  the  method  developed  herein. 

9.  VISCO  results  using  the  Bodner-Partom  material  model  were  very 
similar  to  the  VISCO  results  incorporating  the  Malvern-Norton 
superposition  model  for  the  center  cracked  plate,  especially  for 
times  greater  than  200  seconds  after  load  application.  Also  for 
times  greater  than  200  -.econds,  Norton's  law  alone  in  VISCO  was 
very  similar  to  VISCO  results  using  the  Bodner-Parton  material 
model . 

10.  Displacements  and  the  associated  crack  growth  were  found  to  be 
significantly  dependent  on  prior  deformation  history.  Prior 
deformation  history  became  very  important  in  the  case  of  Test  8 
where  no  fatigue  precracking  was  done  between  creep  crack  growth 
tests. 

The  above  advancements  in  the  understanding  of  creep  crack  growth 
behavior  at  elevated  temperature  are  especially  suited  for  aiding  future 
slow  crack  growth  tests  for  determining  the  threshold  load  levels  for 
creep  crack  growth.  In  addition  the  crack  growth  criteria  investigations 
provide  significant  progress  towards  life  prediction  of  actual  turbine 
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The  following  additional  research  work  is  recommended  to  further 
the  life  prediction  capability  developed  in  the  present  work. 

1.  Use  present  approach  to  analyze  other  test  specimen  geometries 
(e.g.,  compact  tension  specimen)  to  determine  dependence/ Independence 

of  results  on  specimen  geometry  and  also  assess  repeatability  and  material 
data  scatter. 

2.  As  computers  become  faster  and  more  efficient  finite  element 
techniques  are  developed,  a  three-dimensional  analysis  of  creep  crack 
growth  should  be  accomplished  to  correctly  model  through  the  thickness 
variations  in  test  specimen  behavior.  As  a  first  step  in  this  direction, 
plane  strain  analyses  similar  to  the  present  work  might  provide  addi¬ 
tional  insight  into  creep  crack  growth  behavior. 

3.  Future  work  needs  to  include  cyclic  or  engine  spectrum  loading 
conditions. 

4.  Additional  material  characterization  test  data  is  needed  in 
general  for  IN-100.  Tests  providing  this  data  should  be  done  at  several 
temperatures  such  that  the  constitutive  model  could  be  further  developed 
and  include  temperature  dependence. 

5.  Environmental  effects  should  be  further  researched  by  perform¬ 
ing  creep  crack  growth  testing  in  several  different  environments  such 
as  vacuum,  inert,  salt  spray,  high  sulfur  content,  variable  oxygen 
partial  pressures,  etc. 

In  summary  the  technique  developed  herein  is  worth  exploring  further 
in  that  it  has  potential  for  providing  information  on  crack  growth  rate 
behavior  in  engine  materials  under  typical  operating  conditions  which 
might  not  be  readily  obtained  using  conventional  techniques.  This 
Information  in  turn,  is  necessary  in  order  to  implement  a  retirement 
for  cause  philosophy  for  U.S.  Air  Force  jet  engine  components. 
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APPENDIX  A 

FINITE  ELEMENT  FORMULATION 

In  the  following  sections  the  basic  concepts  and  equations  used  in 
the  finite-element  analysis  of  elastic  and  elastic-plastic  materials  are 
briefly  reviewed.  These  equations  provide  background  for  the  elastic- 
viscoplastic  nonlinear  finite-element  computer  program  development. 

The  basic  philosophy  of  the  finite-element  method  (Reference  31)  is 
that  an  approximate  solution  to  a  complicated  problem  can  be  obtained  by 
subdividing  the  region  of  interest  into  a  finite  number  of  discrete  ele¬ 
ments  and  then  choosing  appropriate  relatively  simple  functions  to 
represent  the  solution  within  each  element.  These  functions  are  simple 
compared  to  the  so-called  "exact"  solutions  which  account  for  the  entire 
region  of  interest.  In  this  section  the  equations  associated  with 
representing  a  two-dimensional  body  as  a  finite  number  of  elements  are 
presented.  The  displacements  in  each  element  were  expressed  as  a  simple 
polynomial  and  the  equations  relating  displacements  to  applied  loading 
for  both  plane-stress  and  plane-strain  conditions  are  given. 

1.  DISPLACEMENT  MODEL 

The  displacement  function  used  in  the  displacement  formulation  is 
generally  selected  as  a  polynomial.  The  polynomial  expression  allows 
for  simple  differentiation  and  integration.  Also,  as  the  element  size 
becomes  small,  the  polynomial  expression  permits  a  simple  approximation 
to  the  exact  solution.  A  polynomial  of  infinite  order  corresponds  to 
an  exact  solution.  However,  for  practical  purposes  the  polynomial  must 
be  truncated  to  a  finite  number  of  terms.  Thus,  the  number  of  elements 
in  a  structure  must  be  large  enough  so  that  the  displacement  function 
for  each  element  closely  approximates  the  exact  displacements  in  that 
particular  region. 

In  any  numerical  method,  the  solution  should  converge  to  the  exact 
solution  as  the  size  of  the  elements  become  small.  For  the  displacement 
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formulation,  it  has  been  shown  that  under  certain  conditions  the  solution 
provides  a  lower  bound  to  the  exact  displacements  (Reference  31).  To 
assure  this  convergence  certain  conditions  must  be  satisfied.  First,  the 
displacement  function  must  be  chosen  so  that  rigid  body  displacements  do 
not  cause  straining  of  the  element.  Second,  the  function  must  also  be 
chosen  so  that  a  constant  state  of  strain  is  obtained  as  the  element  size 
approaches  zero.  The  simplest  polynomial  function  which  satisfies  these 
two  requirements  and  also  maintains  displacement  continuity  between  adja¬ 
cent  elements  is  the  linear-displacement  function. 

a.  Displacement  Function 

Figure  A-l  shows  a  typical  triangular  element,  m,  with  nodes  i,  j, 
k  numbered  in  a  counter-clockwise  direction.  The  linear-displacement 
function  which  defines  the  displacements  within  the  element  is  given  by 


U  =  ai  +  a2  X  +  <*3  y 
V  =  a4  +  a5  X  +  a6  y 


(A-l) 


where  the  constants  a.,  are  determined  from  the  six  nodal  displacements 
and  nodal  coordinates  as 
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(A-3) 


where  A  Is  the  In-plane  area  of  the  element.  The  coefficients  a,,  b. , 
m  11 

and  c^  are  given  by 

at  '  XA  -  Vj  (A-4> 
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Figure  A-1.  Constant  Strain  Triangular  Element 
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b1  ■  yj  -  yk  (A-4) 

C1  =  Xk  -  xj 

where  x  and  y  are  coordinates  of  the  nodal  points.  The  other  coefficients 
for  subscripts  "j"  and  "k"  are  obtained  by  cyclic  permutation  of  the 
subscripts  i,  j,  and  k. 


b.  Element  Strain 

The  total  strains  at  any  point  within  an  element  are  defined  in  terms 
of  the  displacement  derivatives  as 
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(A- 5) 


From  Equations  A-l  to  A-5,  the  total  strains  are  written  in  terms  of 
nodal  displacements  and  coordinates  as 


{t)  =  [B]  {U} 


(A— 6) 


where  {U>  is  the  generalized  nodal  displacement  (U}T 
and 


'ViYAV 


The  superscript  T  denotes  the  matrix  transpose. 


(A-7 ) 


2.  ELASTIC  ANALYSIS 


For  linear  elastic  and  isotropic  materials,  the  relationship  between 
stresses  {ol,  strains  (el.  Initial  stresses  {o0>  and  any  initial  strains 
{eQ}  is  given  by 


(a)  * 


[0]  {{el  -  (e0}}  +  {aQ> 


(A-8) 


o 


xy 
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where  [D]  is  the  elastic  material  property  matrix.  The  matrix  [D]  for 
plane-stress  conditions  where  oz  =  axz  =  ayz  «  0  is  given  by 


(A-  9 ) 


where  E  and  v  are  the  modulus  of  elasticity  and  Poisson's  ratio, 
respectively.  For  plane-strain  conditions  where  cz  =  0,  the  elastic 
material  property  matrix  is  given  by 


[D] 


=  EO-y) 

(l+v)(l-2v) 


1  v/l-v 

v^l-v  1 
0  0 


0 

0 


1  -2v 

2(1 -v ) 


(A-10) 


Under  plane-strain  conditions  a  normal  stress  also  exists  and  is  given  by 


a  =  v  (o  +  a  ) 
z  x  y 


( A— 1 1 ) 


a.  Method  of  Solution 

The  equation  which  governs  the  elastic  response  of  a  discretized 
structure  can  be  derived  from  the  principle  of  virtual  work  (Reference 
31)  and  is  given  by 

[K]  {U}  =  { P>  +  {Q}  (A-l 2 ) 

where  [K]  is  the  elastic  stiffness  matrix  of  the  structure,  { U>  is  the 
generalized  displacement  vector,  { P)  is  the  external  applied  load  vector, 
and  (Q)  is  the  force  vector  due  to  the  presence  of  initial  stress  and/or 
initial  strain. 
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The  coefficients  of  the  elastic  stiffness  matrix  are  obtained  from 

[K]  =  E  f [B]T[D]  [B]  dvol  (A- 13) 

m=l  J 

where  the  integration  is  taken  over  the  volume  of  each  element  and  the 
summation  is  over  all  elements  in  the  structure.  The  nodal  forces  due  to 
initial  stresses  are  given  by 

(Q }CTo  =  Y  fl B]T{%>  d  vol  (A-14) 

m=  1 

and  the  nodal  forces  due  to  initial  strains  are  given  by 

M 

{Q}eo  =  Y  /[B]T  CD]  {£o}  d  vo1  (A-15) 

m=l  J 

3.  ELASTIC-PLASTIC  ANALYSIS 

The  application  of  the  finite  element  method  to  problems  involving 
materials  that  obey  linear  constitutive  laws  is  straightforward  because 
the  material  properties  are  constant.  Therefore  only  one  solution  is 
required  to  obtain  displacements  for  the  elastic  structure.  However,  for 
elastic-plastic  problems  the  coefficients  in  the  stiffness  matrix  vary  as 
a  function  of  loading.  Thus,  the  elastic-plastic  displacements  are 
usually  obtained  by  applying  small  load  increments  to  the  structure  and 
updating  the  coefficients  of  the  stiffness  matrix.  Another  technique 
called  the  "residual  force"  method  (Reference  11)  avoids  modifications 
of  the  stiffness  matrix  by  adding  on  a  so-called  plastic  load  vector 
to  the  force  side  of  the  equilibrium  equation  (i.e.,  Equation  A-12). 

Only  the  residual  force  method  will  be  discussed  herein. 

a.  Yield  Criterion 

In  any  elastic-plastic  material  the  elastic  formulation  can  be  used 
prior  to  plastic  yielding.  Thereafter  it  is  necessary  to  have  a  yield 
criterion  to  determine  the  state  of  stress  at  which  yielding  occurs.  The 
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von  Mises  yield  criterion  is  one  of  the  most  widely  used.  This  criterion 
assumes  yielding  is  caused  by  the  maximum  distortion  energy  (Reference 
64).  The  yield  criterion  for  plane  stress  conditions  is  given  by 

F  =  F{a)  =  [ax2  +  cry2  -  +  3ox2j^-  a  (A-16) 

and  for  plane  strain 

F  =  F{o)  =  [ax2  +  ay2  +  az2  -  -  o)(o1  -  cyjz  +  3crx2]  -  a  (A-17) 

where  a  is  the  uniaxial  yield  stress.  If  the  state  of  stress  is  such 
that  F  <  0,  the  material  is  still  in  the  elastic  range.  When  F  =  0,  a 
plastic  state  is  obtained  and  one  of  the  flow  theories  of  plasticity 
must  be  employed  to  determine  subsequent  plastic  behavior  under  increas¬ 
ing  stress  or  strain. 

b.  Flow  Theory 

One  of  the  basic  assumptions  in  the  theory  of  plasticity  is  that 
the  total  strain  (e)  or  total  strain  increment  {del  can  be  decomposed 
into  elastic  and  plastic  strain  components  as  follows: 

{ e }  =  {cE}  +  (cp>  (A- 18) 

or,  incrementally, 

(de)  =  (deE)  +  {deP}  (A-19) 

In  the  incremental  theory  of  plasticity  the  plastic  strain  increment 
vector  {deP}  is  a  function  of  the  current  state  of  stress  and  is  related 
to  the  yield  criterion  through  the  Theory  of  Plastic  Potential  (Refer¬ 
ence  80) 

(dEP)  -  *{s&}  (A-20) 

where  x  is  a  positive  scalar  quantity.  This  flow  law  is  also  written 
in  terms  of  strain  rate 
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In  this  case  x'  has  the  significance  of  the  coefficient  of  viscosity. 
Equations  A-20  and  A-21  are  also  known  as  Drucker's  Normality  Principle 
(Reference  64)  which  by  its  name  specifies  that  the  plastic-strain 
increment  vector  is  to  be  aligned  normal  to  the  yield  surface  in  nine¬ 
dimensional  stress  space.  When  the  von  Mises  yield  criterion  is  used 
with  Equation  A-20  the  resulting  expression  for  {de*3>  is  identical  to 
that  proposed  by  Prandtl  and  Reuss  (Reference  64).  The  total  strain 
increment  vector  can  now  be  written  as 

{del  =  [0]{do)  +  (A-22) 

where  the  elastic  strain  increment  vector  has  been  related  to  the  stress 
increments  {da}  through  the  elasticity  matrix.  Therefore,  if  x  was  known, 
then  the  desired  stress-strain  relation  for  an  elastic-plastic  material 
would  be  obtained.  When  yielding  is  occurring,  the  total  differential 
of  Equation  A-16  or  Equation  A-17  gives 

dF  =  {ff^}}T  {do}  "  d°  =  0  (A-23) 

The  increment  in  yield  stress  da  is  obtained  from  a  uniaxial  tensile 
test  as 

dF  =  dou  =  deP  =  H'  dgP  (A-24) 

where  H1  is  the  slope  of  the  stress-plastic  strain  curve,  dou  is  the 
uniaxial  stress  increment,  d eP  is  the  uniaxial  plastic  strain  increment. 
Using  Equation  A-20  for  the  uniaxial  case  gives  dep  =  x.  Thus  Equation 
A-23  becomes 

{  Ifcy}  {do}  -  H'  x  =  0  (A-25) 
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Eliminating  a  from  Equations  A-22  and  A-25  results  in  an  explicit  expres¬ 
sion  relating  increments  in  stress  to  increments  in  total  strain  (Refer¬ 
ence  44).  This  expression  is 

{do}  =  [0Ep]  {del  (A-26) 


where 

[dep]  =  [0]  m  ["'  +  {!£„,}  m  {!U  ]'’  <A-27> 

The  matrix  [D^p]  is  the  elastic-plastic  matrix  which  replaces  the 

elasticity  matrix  [D]  in  an  incremental  analysis.  For  an  elastic-perfectly 
plastic  material,  H’  is  set  equal  to  zero.  In  general  the  slope  of  the 
uniaxial  stress-plastic  strain  curve,  H',  varies  with  plastic  strain. 
Therefore,  to  relate  a  multiaxial  plastic  strain  state  to  a  uniaxial 
experimental  stress-plastic  strain  curve,  an  effective  plastic  strain 
is  defined  in  incremental  form  as 

d'e  -  VHj  de?j  (A-28> 

4.  ELASTIC-PLASTIC  SOLUTION  TECHNIQUES 

The  procedures  used  to  solve  small  displacement  elastic-plastic 
problems  incrementally  within  a  finite  element  computer  program  may  be 
divided  into  two  categories.  In  one  the  effects  of  plasticity  are 
accounted  for  directly  in  the  stiffness  matrix.  The  second  category 
treats  plastic  behavior  as  an  additional  plastic  load  that  is  combined 
with  applied  or  external  loads  in  the  equilibrium  equation  (i.e.  Equa¬ 
tion  A-12).  These  two  procedures  are  referred  to  as  the  "tangent  modu¬ 
lus"  and  "residual  force"  methods  respectively.  Only  the  residual  force 
method  in  the  form  of  "initial  stress"  and  "initial  strain"  will  be 
summarized  herein. 
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a.  Initial  Stress  Method  (Reference  44) 

The  equation  which  governs  the  response  of  a  discretized  structure 
under  loads  which  cause  plastic  deformation  (Reference  31)  is 

[K]  (Ulj  =  {P}1  +  (A-29) 

where  [K^]  is  the  elastic  stiffness  matrix,  {U}  is  the  generalized  dis¬ 
placement  vector,  { P>  is  the  applied  load  vector,  and  { Q}  is  the  "effec¬ 
tive"  plastic-load  vector  which  accounts  for  elements  in  a  plastic  state. 
The  initial  stress  method  approaches  the  solution  to  an  elastic  plastic 
problem  by  applying  a  series  of  small  load  increments  to  the  structure 
until  the  desired  load  is  reached  (  {P}1  =  { P> 1 ~ ^  +  {dP}).  The  super¬ 
script  i  denotes  the  current  increment  and  i-1  denotes  the  preceding 
increment.  After  each  load  increment  an  iterative  process  is  required 
to  stabilize  the  plastic-load  vector.  The  subscript  I  denotes  the 
current  iteration  and  1-1  denotes  the  preceding  iteration.  During  the 
ith  increment  a  purely  elastic  problem  is  solved  and  the  increments  in 
total  strain  (de)  and  corresponding  elastic  stress  {doE>  are  computed 
from  the  displacement  increments  { dUJ  for  every  element.  Because  of 
the  material  nonlinearity  the  stress  increments  are  not,  in  general 
correct  or  if  the  correct  stress  increment  for  the  corresponding  strain 
increment  is  {da},  then  a  set  of  body  forces  or  plastic-load  vectors 
{ dQ>  caused  by  the  "initial"  stress  {daQ}  =  {da^}  -  {da})  is  required  to 
maintain  the  stress  components  on  the  yield  surface  or  compatible  with 
the  uniaxial  stress-strain  curve.  The  correct  stress  increment  is  com¬ 
puted  with  Equation  A-26.  The  plastic  load  increments  are  computed  from 
Equation  A-14 

M  f  T 

{dQ}  =  El  CB]  {da  }  d  vol  (A-30) 

j  o 

m=l 

Elements  are  in  the  elastic  state  when  {doQ}  =  0.  The  total  plastic-load 
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vector  is  then  computed  as 

{Q>j  =  (Q> j”]  +  (dQ}j  (A-31) 

At  the  second  stage  of  computation  the  new  force  system  (Q}J  is  added  to 

the  applied  load  vector  and  a  new  set  of  displacements  is  obtained. 

Again,  some  of  the  stresses  are  likely  to  exceed  the  yield  criterion 
and  a  new  set  of  plastic-load  increments  are  computed.  This  iteration 
process  is  repeated  until  the  change  in  the  plastic-load  vector  is 
sufficiently  small.  See  Figure  A-2  for  a  uniaxial  schematic  of  this 
iterative  procedure  and  Figure  A-3  for  the  mathematical  algorithm. 

Consider  points  A,  B,  and  C  in  Figure  A-2.  Point  A  is  the  state  prior 
to  the  load  increment.  Point  B  is  the  state  after  the  load  increment 
has  been  applied  and  one  "initial"  stress  iteration  has  been  accomplished. 
Point  C  is  the  state  of  stress  and  strain  sought  after  which  satisfies 
equilibrium  with  external  loads  and  compatibility  with  the  material's 
stress-strain  curve.  Notice  point  B  satisfies  compatibility  but  not 
equilibrium  since  {do}j  <_  (do^lj. 

5.  INITIAL  STRAIN  METHOD 

The  initial  strain  method  parallels  the  initial  stress  method 
somewhat  and  accordingly  this  development  will  begin  just  after  Step 
4  in  the  "Initial  Stress"  algorithm  in  Figure  A-3. 

The  elastic-plastic  material  matrix  [D^p]  is  used  as  follows 

{dePjj  =  {de}j  -[D]"1  [DEp]  {dclj  (A-32) 

This  plastic  strain  increment  {dep} j  is  then  used  to  calculate  a  plastic 
force  vector  increment 

M 

{dQ>j  =  £  /[ B]T  [D]  LdePjj  d  vol  (A-33) 

m=l 
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Figure  A-2.  Uniaxial  "Initial"  Stress  Iteration  Schematic 
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For  Each  Load  Increment  dP  in  Plastic  Range  44 


r  -1 

(  { dP)  t  {dQ>r  -j 

1. 

{dUlT 

*  m 

2. 

(del  j 

=  [B]  {dUlj 

3. 

{dog.} 

r  =  [0]  {de>j 

4. 

(a'lj 

=  to}N1 

+  {do^lj 

5. 

{do)I 

[oEP]i 

{dtlj 

6. 

tdco} 

=  {doE}j 

-  (do)  j 

7. 

{a}  j 

=  {a,l1  - 

(do  } 

o  I 

{  £  1 1 

=  {e}I-l 

+  {delj 

8. 

{dQlj 

fl  b]t 

{do  )T  dvol 

0  1 

9. 

Continue  steps 

1+8  until  {dQ}j 

Figure 

■  A-3. 

Initial 

Stress  Algorithm 
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This  plastic  force  vector  increment  is  added  to  the  external  force  vector 
increment  {dP}  for  the  augmented  global  force  vector  used  in  the  next 
iteration  as  follows: 


{dU)I+1  =  [K]"1  {  (dP)  +  {dQ}j} 

{de)I+1  =  [B]  {dUlI+1 

The  stress  increment  (dal  is  calculated  as  follows 

{do} I+1  =  [D]  {  (de)  -  {dePjj} 

Steps  (1)  -  (4)  of  the  Initial  Stress  algorithm  and  Equations  A- 32  through 
A-36  above  are  repeated  until  compatibility  with  the  materials  stress- 
strain  curve  is  established.  Compatibility  is  shown  to  be  achieved  after 
"n"  iterations  in  Figure  A-4.  Also,  compatibility  would  display  itself 
by  little  or  no  change  in  the  plastic  strain  increment  between  iterations. 
Note  that  equilibrium  is  continually  satisfied  in  this  initial  strain 
method.  This  version  of  the  initial  strain  method  differs  from  Marcal 
(Reference  45)  by  the  fact  that  iterations  within  a  load  increment  are 
done  globally  rather  than  within  each  element  as  the  "constant  strain" 
method  of  iteration  implies  in  Marcal 's  paper. 


(A— 34) 
(A-35) 

(A-36) 
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Figure  A-4.  Uniaxial  Initial  Strain  Iteration  Schematic 
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APPENDIX  B 

ITERATIVE  SOLUTION  TECHNIQUE 
FOR  NODE  POINT  DISPLACEMENTS 

The  matrix  equation  which  governs  the  response  of  a  discretized 
structure  is 

[K]  (U)  =  (P)  (B-l ) 

where  [K]  is  a  symmetric  positive-definite  n  x  n  matrix,  (U>  is  the  unknown 
node  point  displacement  vector,  and  (P)  is  a  known  load  vector.  In  the 
finite  element  method  for  structural  analysis,  the  matrix  [K]  is  usually 
highly  banded  and  if  stored  in  compacted  form  (i.e.,  only  nonzero  terms 
retained)  requires  much  less  space  in  the  computer  than  the  product  n  x  n 
reflects.  Also,  if  there  are  changing  boundary  conditions,  such  as  free¬ 
ing  nodes  to  simulate  crack  growth,  then  the  [K]  matrix  must  be  recomputed. 
A  solution  technique  that  works  well  with  compacted  [K]  matrices  and 
conveniently  admits  boundary  condition  changes,  is  the  Gauss-Seidel 
iterative  technique  with  over-relaxation  (Reference  14).  This  technique 
may  be  implemented  in  the  following  manner  (Reference  87).  Consider 
Equation  B-l  rewritten  as 


[K]  {U>  = 


S  X  X  !  S  X  Y 


S  Y  X  ;  S  Y  Yj 


(B-2) 


where  U  and  U  represent  node  point  displacement  vectors  in  the  x  and  y 
x  y 

direction  respectively,  P  and  P  represent  the  node  force  vectors  in  the 

x  y 

x  and  y  direction  respectively.  The  submatrices  SXX,  SXY,  SYX,  and  SYY 
in  the  matrix  [K]  have  dimensions  ^  x  but  due  to  their  bandedness 
can  be  compacted  to  a  matrix  which  is  ^  x  9.  The  dimension  9  minus  1 
reflects  how  many  adjacent  nodes  can  be  connected  to  any  given  node. 
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This  is  not  very  restrictive  since  triangular  finite  elements  develop 
undesirable  aspect  ratios  if  there  are  more  than  eight  nodes  connected 
to  any  one  node. 

Appropriate  terms  of  the  matrix  [K]  are  retrieved  from  the  compacted 
submatrices  with  the  help  of  a  matrix  NP  and  a  vector  NAP.  The  vector 

NAP  has  the  dimension  j  and  the  matrix  NP  has  the  dimension  ^  x  9.  The 

Ith  component  of  NAP  or  NAP(I)  stores  the  number  of  adjacent  node  points 
connected  to  node  point  I.  The  (I,  J)  component  of  NP  or  NP(I,J)  stores 
the  address  of  the  terms  in  the  submatrices  associated  with  the  Jth 
adjacent  node  point  connected  to  node  point  I.  Note  that  for  node  I,  J 
may  go  from  1  to  NAP(I). 


Consider  the  governing  equation  for  node  point  I  displacements 
written  as 


SXX  (I,  1)  SXY  (I,  1)  (  Ux  (I) 

SYX  (I,  1)  SYY  (I,  1)  Uy  (I) 


P  (I) 
x  '  ' 

Py  (I) 


SXX  (I,J) 
SYX  (I,J) 


J=2 


SXY  (I,J)1  f  Ux  (J)) 

SYY  ( I , J )J  (  Uy  (J)J  (B-3) 


If  the  right-hand  side  of  Equation  B-3  is  defined  as  the  vector 
then  solving  for  the  displacements  at  node  I  yields 


ux  (I) 

uy  (i) 


SXX  (I,  1) 
SYX  (I,  1) 


SXY  (I,  1) 
SYY  (I,  1) 


(B-4) 


Note  that  the  matrix  to  be  inverted  in  Equation  B-4  is  only  a  2  x  2. 
Also  since  this  is  the  only  place  these  terms  of  the  submatrices  are 
used  this  2x2  may  be  inverted  and  its  components  stored  in  their 
original  submatrix  locations. 
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To  incorporate  an  over-relaxation  factor.  Equation  B-4  is  modified 
as  follows 


i 


AU. 


m 

P  - 

-1, 

m 

(I) 

SXX  (I,  1)  SXY  (I,  1) 

FRX 

.  jv1’ 

(I) 

SYX  (I,  1)  SYY  (I,  1) 

|  FRY 

1 u,  <» 

(B-5) 


where  the  superscripts  m  and  m-1  refer  to  iteration  number.  The  left- 
hand  side  of  Equation  B-5  is  the  change  in  displacements  between  iterations 
without  applying  an  over-relaxation  factor.  But  the  new  total  displace¬ 
ments  for  iteration  m  using  an  over-relaxation  factor  are 


Ux  (I))m  (Ux(I))m_1  AUx  (I) 

Uy  (I)j  =  |uy  (I)  j  +  X  FAC  AUy  (I) 


(B-6) 


where  XFAC  is  the  over-relaxation  factor  which  normally  ranges  from  1.8 
to  1.9  for  structural  analysis. 


Convergence  of  these  iterations  is  checked  by  computing  an  effective 
force  unbalance  term,  SUM,  defined  as 

n/2 

SUM  =  £  [|AU  (I)  *  SXX  (I,l)|  +  |  AU  (I)  *  SYY  (1 ,1 )  |  ]  (B-7) 

1  =  1 

If  SUM  becomes  less  than  a  specified  small  value,  e,  iterations  are 
stopped  and  the  node  point  displacement  solution  is  obtained.  The  value 
of  e  is  chosen  based  on  examining  solutions  for  various  sizes  of  e. 

A  good  starting  value  for  e  is  one  tenth  of  the  applied  load.  The  final 
e  is  then  chosen  based  on  the  amount  of  accuracy  desired. 

Displacement  boundary  conditions  are  easily  input  to  this  solution 
routine  by  simply  specifying  the  desired  node  point  displacements  and  then 
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require  the  iteration  algorithm  to  skip  over  the  node  displacement  equa¬ 
tions  which  have  fixed  displacements.  Likewise,  if  a  fixed  node  is 
released  during  the  analysis  such  as  for  modeling  crack  growth,  this 
node's  displacement  equation  may  be  reactivated  in  the  iteration 
algorithm. 

Convergence  of  this  solution  technique  (i.e.,  when  SUM  <_e)  is 
dependent  on  the  initial  guess  for  the  node  displacements.  Usually  for 
convenience  all  unknown  displacements  are  initialized  at  zero.  However, 
for  each  succeeding  solution,  such  as  in  a  nonlinear  incremental  analysis, 
much  better  initial  displacement  values  are  available  from  the  prior  solu¬ 
tion.  These  initial  displacements  from  prior  solutions  significantly 
reduce  the  number  of  iterations  to  convergence. 
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APPENDIX  C 

THE  J  AND  C*  INTEGRALS 

1.  J  INTEGRAL 

Rice's  0  integral  is  defined  as 

J  =  /  [  W  (e.p  dy  -  Ti  f^-ds  ]  (C-l) 

where 

/Emn 

„  °u  deij  (c‘z) 

and  r  is  a  closed  loop  around  the  crack  tip  as  shown  in  Figure  2. 
Expanding  J  and  integrating  along  a  rectilinear  path  in  the  x  and  y 
directions  results  in 

J  =  /  C  (W  -°x  «  '  °xy  «  }  dy  +  (oxy  U  +  °y  M  }  dx  ]  (C'3) 

r 

The  following  describes  a  numerical  procedure  for  calculating  the  J 
integral  with  a  finite  element  program  that  incorporates  constant  strain 
triangles.  Consider  Figure  C-l  which  is  a  region  of  elements  taken  from 
a  finite  element  model  of  a  cracked  geometry.  Paths  1  and  2  are  two 
possible  paths.  The  contribution  to  Path  1  from  Element  2  is 

4?™  2  ■  tW<2)  -  °x<2>  |f  (2)  -  °xy  (2)  |f  (2)  ]  «  (C-4) 

where  the  number  inside  parenthesis  refers  to  the  element  number  these 
values  came  from.  A  similar  contribution  would  come  from  Element  5  for 
Path  1.  Path  2  which  also  runs  through  Element  2  would  also  have  the 
above  contribution  from  Element  2  but  as  Path  2  turns  and  runs  along 
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y 

v 


Figure  C- 


1.  J  Integral  Paths  within  a  Constant  Strain  Triangul 
Finite  Element  Model 
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the  border  of  Elements  6  and  1  the  following  is  the  contribution 


Aj 


Path  2 

Elements 
1  &  6 


'  t»xy  <6>  S  <6>  *  °y  («  £ 

*  [°xy  <’>  £  (1)  +  °y  (’>£(')]  (T) 


(C-5) 


Notice  in  this  case  an  average  of  the  stresses  and  strains  in  the  two 
elements  is  taken  by  effectively  running  half  the  element  length  in 
Element  6  and  half  in  Element  1. 


The  strain  energy  term  in  Equation  C-4  is  calculated  as 
follows 


W(2)  = 


1 

2  °ij 


(C-6) 


Although  these  equations  are  shown  for  specific  element  numbers  their 
form  is  used  for  all  elements  along  the  J  integral  path. 


2.  C*  INTEGRAL 

The  C*  Integral  (Reference  23)  defined  as 


where 


C* 


§  [W*  dy  -  T.  ds] 

r 


W* 


(C-7) 


(C-8) 


and  r  is  the  same  type  of  path  as  for  the  J  integral.  C*  may  be  obtained 
by  replacing  strain  and  displacement  in  the  J  integral  with  strain  rate 
and  displacement  rate  respectively.  The  rational  behind  this  is  based 
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on  the  assumption  that  the  material  being  analyzed  behaves  as  a  creeping 
solid  such  that 

ep  =  y  (a)e  (C-9) 

or  in  a  multiaxial  sense 


where  W*  was  defined  in  Equation  C-8.  Notice  Equation  C-9,  which  is 
considered  a  plastic  strain  rate,  makes  no  provision  for  admitting  elastic 
strain  rates.  Hence  the  claim  of  path  independence  for  the  C*  integral 
can  only  be  approached  in  a  realistic  elastic-plastic  material  once  the 
stresses  have  reached  a  steady  state  value  (i.e.,  a.,  =  0).  With  this 

restriction  on  stresses,  C*  may  also  be  defined  as  j  or  the  time  rate  of 
change  of  the  J  integral.  In  general  j  would  include  the  time  derivatives 
of  stress  and  traction  but  if  they  are  restricted  to  zero  then  j  is  equal 
to  C*  for  a  creeping  solid.  Also  with  o..  equal  to  zero.  Equation  C-8 

•  J 

may  be  directly  integrated  to 

“*  ■  °ij  <c-"> 

Using  this  form  of  W*,  C*  may  be  numerically  calculated  in  the  same 
fashion  as  the  J  integral  with  the  exception  of  using  strain  rates  and 
displacement. 
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APPENDIX  D 

DETERMINATION  OF  BODNER 
MATERIAL  MODEL  CONSTANTS 

This  appendix  describes  how  each  of  the  material  constants  for  the 
Bodner-Partom  material  model  can  be  determined.  Normally  high  costs 
and  material  shortages  prevent  obtaining  more  than  one  or  two  stress- 
strain  curves  and  the  same  number  of  creep  tests.  Ideally,  to  best 
characterize  the  material,  several  stress-strain  curves  should  be 
generated  over  a  wide  range  of  strain  rates,  similar  to  Figure  5. 
Likewise,  several  creep  tests  should  be  performed  over  a  wide  spectrum 
of  stress  levels. 


In  general  the  Bodner  material  constants  are  dependent  on  tempera¬ 
ture.  However,  the  temperature  dependence  is  suppressed  by  performing 
the  material  characterization  tests  ( i . e. ,  stress-strain  and  creep)  at 
the  same  temperature  that  the  Bodner  model  will  be  applied. 


To  determine  the  Bodner  constants  from  uniaxial  test  data,  the 
Bodner  equations  are  written  in  uniaxial  form  as  follows:  (the  total 
strain  rate  is  the  sum  of  plastic  and  elastic  strain  rates) 


•  _  o  ,  ‘P 
£'E  +  e 


and  the  plastic  strain  rate  for  uniaxial  tension  is 


where 


eP  =k  Do exp[-  "r ] 


Z  =  Z1  -  (Z1  -  Z0)  exp[-mWp] 


V*£P  +  if^z) 


(D-l) 


(D-2) 

(D-3) 

(D-4) 


168 


AFWAL-TR-80-4140 


Z 


rec 


(D-5) 


.3 

over  small  time  intervals,  such  as  a  typical  stress-strain  test  at  10 
sec’1  strain  rate,  the  recovery  term,  Zrec,  may  be  ignored  and  then 

o  e  >  (D-6) 

P 


Since  test  data  can  be  resolved  into  the  forms,  ep  and  0,  Equation 
D-2  is  solved  for  Z  which  then  is  a  function  of  ep  and  a  as  follows 


Z  = 


(D-7) 


The  viscoplastic  material  constants  in  these  equations  are  broken 
into  "short  time  response"  and  "creep"  groupings  for  determination. 


l.  SHORT  TIME  RESPONSE  CONSTANTS 

The  short  time  response  constants  for  the  Bodner  model  are  DQ ,  n, 

m,  ZQ,  ly  These  constants  are  primarily  determined  by  using  stress- 

strain  test  data. 

4  -1 

The  constant  0Q  is  normally  assigned  the  value  of  10  sec  .  For 
high  strain  rate  applications,  DQ  may  be  set  higher  (e.g.,  106  sec’1) 
which  would  result  in  small  changes  to  the  other  constants. 

The  constant  n  is  directly  related  to  the  model's  strain  rate 
sensitivity.  High  n  values  reflect  low  strain  rate  sensitivity  and 
vice  versa.  Changes  in  n  affect  the  stress  values  for  a  given  strain 
rate  by  shifting  up  or  down  the  family  of  stress-strain  curves,  but  the 
shape  of  the  curve  is  preserved.  The  value  of  n  is  determined  in  an 
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iterative  fashion.  The  first  estimate  of  n  should  be  between  1  and  10 
based  on  past  material  modeling  efforts  (Reference  8).  Plots  of  Z  versus 
Wp  are  then  made  by  simultaneous  use  of  each  stress-strain  curve.  Equa¬ 
tions  D-l,  D-6,  and  D-7.  The  value  of  n  is  then  adjusted  with  the  objec¬ 
tive  of  making  all  Z  versus  Wp  curves  from  each  stress-strain  test  fall 

on  top  one  another.  This  value  of  n  will  then  satisfy  the  requirement 
that  Z  is  a  single  value  function  of  Wp  as  given  in  Equation  D-3. 

The  first  approximation  of  ZQ  can  also  be  determined  from  Equation 

D-7.  A  small  value  of  ep  (e.g.,  10-6  sec ""*)  and  the  lowest  apparent 

initial  yield  stress  from  the  stress-strain  test  data  are  substituted 
into  Equation  D-7.  The  resulting  value  of  Z  is  defined  as  ZQ.  Note 

that  ZQ  is  the  primary  constant  that  determines  the  stress  level  at  which 
significant  plastic  straining  (i.e.,  eP  ^  10“6  sec-1)  begins. 

The  constants  Z-j  and  m  are  determined  by  rewriting  Equation  D-3  as 

In  (Z]  -  Z)  =  In  (Z1  -  ZQ)  -  m  Wp  (D-8) 

An  iterative  process  is  now  begun  to  determine  Z^.  The  first  estimate 

of  Z,  should  be  larger  than  Z  (e.g.,  1.5  Z  )  since  Z.  is  the  maximum 

value  for  Z.  By  incorporating  this  estimate  of  Z-| ,  a  plot  of  In  (Z^-Z) 

versus  Wp  is  made  based  again  on  stress-strain  test  data.  Equations  D-l, 

D-6,  and  D-7.  This  should  approximate  a  straight  line  whose  slope  is 
the  constant  m  and  extrapolation  of  the  line  to  Wp  ~  0  provides  a  value 

of  ZQ.  If  this  ZQ  obtained  graphically  does  not  agree  with  the  previous 

value  for  ZQ,  adjust  Z-j  accordingly  and  reiterate. 

These  values  for  n,  m,  ZQ,  and  Z-|  which  primarily  govern  the 
short-term  stress-strain  behavior  should  be  input  to  a  computer  program 


170 


AFWAL-TR-80-41 40 


that  numerically  integrates  Equations  D-l  through  D-4.  Stress-strain 
predictions  from  the  model  should  be  made  for  each  experimental  strain 
rate  to  see  how  good  the  fit  is.  Only  minor  adjustments  to  the  short- 
time-response  constants  should  be  necessary  to  best  fit  test  data. 


2.  CREEP  CONSTANTS 

The  creep  constants  for  the  Bodner  model  are  Z^,  A,  and  r.  These 

constants  are  determined  based  on  data  from  at  least  two  creep  tests  at 
different  stress  levels. 


During  secondary  creep,  when  plastic  strain  rate  is  approximately 
constant,  it  follows  from  Equation  D-7  that  Z  must  also  be  a  constant. 
In  addition,  if  Z  is  a  constant  it  follows  from  Equation  D-3  that  W 

P 

is  a  constant  which  makes  =  0.  Hence,  with  Wp  =  0  in  Equation  D-4 
and  combining  with  Equation  D-5 


O  £ 


P  = 


-  Z) 


( D— 9 ) 


The  constant  Z^  represents  the  minimum  value  for  Z  corresponding 

to  secondary  creep.  If  creep  occurs  below  the  apparent  yield  stress 
implied  by  ZQ,  the  value  of  Z-j  must  be  less  than  ZQ.  Moreover,  Z. 

must  be  less  than  or  equal  to  the  smallest  value  of  Z  from  Equation  D-7 
when  using  creep  test  data  (i.e.,  and  a).  After  selecting  a  value 
for  Z.,  the  constants  A  and  r  can  be  determined  after  rewriting  Equation 

D-9  in  terms  of  natural  logarithms 


1  n  (a 


e  P)  +  ln[m(Z-|-Z)]  =  ln(A  Z , )  +  r  In 


(D-l 0 ) 


Stress  and  plastic  strain  rate  are  substituted  into  Equation  D-l 0  along 
with  the  appropriate  Z  from  Equation  D-7  and  creep  test  data.  With  data 
from  creep  tests  at  two  different  stress  levels  two  linear  equations  in 
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terms  of  the  two  unknowns  A  and  r  are  developed.  The  constants  A  and  r 
are  then  determined  by  simultaneous  solution  of  these  two  equations 
based  on  Equation  D- 1 0  and  creep  test  data. 

The  complete  Bodner  model  should  then  be  tried  out  in  a  computerized 
numerical  integration  scheme  in  an  effort  to  make  final  adjustments 
to  the  constants  for  a  best  fit  to  the  test  data. 

It  became  apparent  to  the  author  that  in  working  with  this  material 
model  a  sensitivity  study  is  required  in  order  to  see  the  effect  of 
changes  in  the  respective  constants  on  creep  crack  growth.  Yet  it  is 
felt  that  the  results  and  conclusions  presented  in  the  main  body  of  this 
report  provide  realistic  directions  and  trends  of  creep  crack  growth. 


